Biomineralization Toolkit Analysis
Kristen Brown, modified ZD
1/30/2024
Biomineralization Toolkit Analysis
sessionInfo() #provides list of loaded packages and version of R. ## R version 4.3.2 (2023-10-31)
## Platform: aarch64-apple-darwin20 (64-bit)
## Running under: macOS Ventura 13.0
##
## Matrix products: default
## BLAS: /Library/Frameworks/R.framework/Versions/4.3-arm64/Resources/lib/libRblas.0.dylib
## LAPACK: /Library/Frameworks/R.framework/Versions/4.3-arm64/Resources/lib/libRlapack.dylib; LAPACK version 3.11.0
##
## locale:
## [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
##
## time zone: America/New_York
## tzcode source: internal
##
## attached base packages:
## [1] stats graphics grDevices utils datasets methods base
##
## loaded via a namespace (and not attached):
## [1] digest_0.6.34 R6_2.5.1 fastmap_1.1.1 xfun_0.41
## [5] cachem_1.0.8 knitr_1.45 htmltools_0.5.7 rmarkdown_2.25
## [9] lifecycle_1.0.4 cli_3.6.2 sass_0.4.8 jquerylib_0.1.4
## [13] compiler_4.3.2 rstudioapi_0.15.0 tools_4.3.2 evaluate_0.23
## [17] bslib_0.6.1 yaml_2.3.8 rlang_1.1.3 jsonlite_1.8.8First, load the necessary packages.
# load libraries - notes show the install command needed to install (pre installed)
library(goseq)## Loading required package: BiasedUrn## Loading required package: geneLenDataBase## library(dplyr)##
## Attaching package: 'dplyr'## The following objects are masked from 'package:stats':
##
## filter, lag## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, unionlibrary(forcats)
library(ggplot2)
library(gridExtra)##
## Attaching package: 'gridExtra'## The following object is masked from 'package:dplyr':
##
## combinelibrary(tidyr)
library(grDevices)
library(reshape2)##
## Attaching package: 'reshape2'## The following object is masked from 'package:tidyr':
##
## smithslibrary(Rmisc)## Loading required package: lattice## Loading required package: plyr## ------------------------------------------------------------------------------## You have loaded plyr after dplyr - this is likely to cause problems.
## If you need functions from both plyr and dplyr, please load plyr first, then dplyr:
## library(plyr); library(dplyr)## ------------------------------------------------------------------------------##
## Attaching package: 'plyr'## The following objects are masked from 'package:dplyr':
##
## arrange, count, desc, failwith, id, mutate, rename, summarise,
## summarizelibrary(ggpubr)##
## Attaching package: 'ggpubr'## The following object is masked from 'package:plyr':
##
## mutatelibrary(tibble)
library(gridExtra)
library(tidyr)
library(zoo)##
## Attaching package: 'zoo'## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numericlibrary(ComplexHeatmap)## Loading required package: grid## ========================================
## ComplexHeatmap version 2.16.0
## Bioconductor page: http://bioconductor.org/packages/ComplexHeatmap/
## Github page: https://github.com/jokergoo/ComplexHeatmap
## Documentation: http://jokergoo.github.io/ComplexHeatmap-reference
##
## If you use it in published research, please cite either one:
## - Gu, Z. Complex Heatmap Visualization. iMeta 2022.
## - Gu, Z. Complex heatmaps reveal patterns and correlations in multidimensional
## genomic data. Bioinformatics 2016.
##
##
## The new InteractiveComplexHeatmap package can directly export static
## complex heatmaps into an interactive Shiny app with zero effort. Have a try!
##
## This message can be suppressed by:
## suppressPackageStartupMessages(library(ComplexHeatmap))
## ========================================library(circlize)## ========================================
## circlize version 0.4.15
## CRAN page: https://cran.r-project.org/package=circlize
## Github page: https://github.com/jokergoo/circlize
## Documentation: https://jokergoo.github.io/circlize_book/book/
##
## If you use it in published research, please cite:
## Gu, Z. circlize implements and enhances circular visualization
## in R. Bioinformatics 2014.
##
## This message can be suppressed by:
## suppressPackageStartupMessages(library(circlize))
## ========================================library(GSEABase)## Loading required package: BiocGenerics##
## Attaching package: 'BiocGenerics'## The following object is masked from 'package:gridExtra':
##
## combine## The following objects are masked from 'package:dplyr':
##
## combine, intersect, setdiff, union## The following objects are masked from 'package:stats':
##
## IQR, mad, sd, var, xtabs## The following objects are masked from 'package:base':
##
## anyDuplicated, aperm, append, as.data.frame, basename, cbind,
## colnames, dirname, do.call, duplicated, eval, evalq, Filter, Find,
## get, grep, grepl, intersect, is.unsorted, lapply, Map, mapply,
## match, mget, order, paste, pmax, pmax.int, pmin, pmin.int,
## Position, rank, rbind, Reduce, rownames, sapply, setdiff, sort,
## table, tapply, union, unique, unsplit, which.max, which.min## Loading required package: Biobase## Welcome to Bioconductor
##
## Vignettes contain introductory material; view with
## 'browseVignettes()'. To cite Bioconductor, see
## 'citation("Biobase")', and for packages 'citation("pkgname")'.## Loading required package: annotate## Loading required package: AnnotationDbi## Loading required package: stats4## Loading required package: IRanges## Loading required package: S4Vectors##
## Attaching package: 'S4Vectors'## The following object is masked from 'package:plyr':
##
## rename## The following object is masked from 'package:tidyr':
##
## expand## The following objects are masked from 'package:dplyr':
##
## first, rename## The following object is masked from 'package:geneLenDataBase':
##
## unfactor## The following object is masked from 'package:utils':
##
## findMatches## The following objects are masked from 'package:base':
##
## expand.grid, I, unname##
## Attaching package: 'IRanges'## The following object is masked from 'package:plyr':
##
## desc## The following objects are masked from 'package:dplyr':
##
## collapse, desc, slice##
## Attaching package: 'AnnotationDbi'## The following object is masked from 'package:dplyr':
##
## select## Loading required package: XML## Loading required package: graph##
## Attaching package: 'graph'## The following object is masked from 'package:XML':
##
## addNode## The following object is masked from 'package:circlize':
##
## degree## The following object is masked from 'package:plyr':
##
## joinlibrary(data.table)##
## Attaching package: 'data.table'## The following object is masked from 'package:IRanges':
##
## shift## The following objects are masked from 'package:S4Vectors':
##
## first, second## The following objects are masked from 'package:reshape2':
##
## dcast, melt## The following objects are masked from 'package:dplyr':
##
## between, first, lastlibrary(stringr)##
## Attaching package: 'stringr'## The following object is masked from 'package:graph':
##
## boundarylibrary(GenomicRanges)## Loading required package: GenomeInfoDblibrary(rtracklayer)
library(rrvgo)Load in data
library(rtracklayer)
gff<-rtracklayer::import("../../data/Pocillopora_acuta_HIv2.genes_fixed.gff3")
gff<-as.data.frame(gff)
dim(gff) # 478988 9## [1] 478988 13names(gff) ## [1] "seqnames" "start" "end" "width"
## [5] "strand" "source" "type" "score"
## [9] "phase" "ID" "transcript_id" "gene_id"
## [13] "Parent"transcripts <- subset(gff, type == "transcript")
transcripts_gr <- makeGRangesFromDataFrame(transcripts, keep.extra.columns=TRUE) #extract length information
transcript_lengths <- width(transcripts_gr) #isolate length of each gene
seqnames<-transcripts_gr$ID #extract list of gene id
lengths<-cbind(seqnames, transcript_lengths)
lengths<-as.data.frame(lengths) #convert to data frame
dim(transcripts) #33730 13## [1] 33730 13kegg <- read.delim("../../data/Pocillopora_acuta_HIv2.genes.KEGG_results.txt",header = FALSE)
kegg <- as.data.frame(kegg)
colnames(kegg)[1] <- "gene_id"
colnames(kegg)[2] <- "KEGG_new"
head(kegg)## gene_id KEGG_new
## 1 Pocillopora_acuta_HIv2___RNAseq.g24143.t1a
## 2 Pocillopora_acuta_HIv2___RNAseq.g24143.t1b K22584
## 3 Pocillopora_acuta_HIv2___RNAseq.g22918.t1
## 4 Pocillopora_acuta_HIv2___RNAseq.g18333.t1 K03386
## 5 Pocillopora_acuta_HIv2___RNAseq.g7985.t1
## 6 Pocillopora_acuta_HIv2___RNAseq.g13511.t1eggnog <- read.delim("../../data/Pocillopora_acuta_HIv2.genes.EggNog_results.txt")#this file contains all of the go terms, descriptions, kegg, etc
eggnog <- plyr::rename(eggnog, c("X.query"="gene_id"))
head(eggnog,2)## gene_id seed_ortholog evalue score
## 1 Pocillopora_acuta_HIv2___RNAseq.g24143.t1a 45351.EDO48725 2.16e-120 364
## 2 Pocillopora_acuta_HIv2___RNAseq.g18333.t1 45351.EDO38694 3.18e-123 355
## eggNOG_OGs
## 1 COG0620@1|root,KOG2263@2759|Eukaryota,38GZS@33154|Opisthokonta,3BNKS@33208|Metazoa
## 2 COG0450@1|root,KOG0852@2759|Eukaryota,38B9P@33154|Opisthokonta,3BGS4@33208|Metazoa
## max_annot_lvl COG_category
## 1 33208|Metazoa E
## 2 33208|Metazoa O
## Description Preferred_name
## 1 Cobalamin-independent synthase, Catalytic domain -
## 2 negative regulation of male germ cell proliferation PRDX4
## GOs
## 1 -
## 2 GO:0000003,GO:0001775,GO:0002252,GO:0002263,GO:0002274,GO:0002275,GO:0002283,GO:0002366,GO:0002376,GO:0002443,GO:0002444,GO:0002446,GO:0003006,GO:0003674,GO:0003824,GO:0004601,GO:0005488,GO:0005515,GO:0005575,GO:0005576,GO:0005615,GO:0005622,GO:0005623,GO:0005737,GO:0005783,GO:0005790,GO:0005829,GO:0006082,GO:0006457,GO:0006464,GO:0006468,GO:0006520,GO:0006575,GO:0006793,GO:0006796,GO:0006807,GO:0006810,GO:0006887,GO:0006915,GO:0006950,GO:0006952,GO:0006955,GO:0006979,GO:0007154,GO:0007165,GO:0007249,GO:0007252,GO:0007275,GO:0007276,GO:0007283,GO:0007548,GO:0008150,GO:0008152,GO:0008219,GO:0008285,GO:0008379,GO:0008406,GO:0008584,GO:0009056,GO:0009266,GO:0009409,GO:0009605,GO:0009607,GO:0009617,GO:0009628,GO:0009636,GO:0009893,GO:0009966,GO:0009967,GO:0009987,GO:0010467,GO:0010604,GO:0010646,GO:0010647,GO:0010941,GO:0010942,GO:0010950,GO:0010952,GO:0012501,GO:0012505,GO:0016043,GO:0016192,GO:0016209,GO:0016310,GO:0016491,GO:0016684,GO:0016999,GO:0017001,GO:0017144,GO:0019222,GO:0019471,GO:0019538,GO:0019725,GO:0019752,GO:0019953,GO:0022414,GO:0022417,GO:0023051,GO:0023052,GO:0023056,GO:0030141,GO:0030162,GO:0030198,GO:0031323,GO:0031325,GO:0031410,GO:0031974,GO:0031982,GO:0031983,GO:0032268,GO:0032270,GO:0032501,GO:0032502,GO:0032504,GO:0032940,GO:0033554,GO:0034774,GO:0035556,GO:0036211,GO:0036230,GO:0042119,GO:0042127,GO:0042221,GO:0042592,GO:0042737,GO:0042742,GO:0042743,GO:0042744,GO:0042802,GO:0042803,GO:0042981,GO:0043062,GO:0043065,GO:0043067,GO:0043068,GO:0043085,GO:0043170,GO:0043207,GO:0043226,GO:0043227,GO:0043229,GO:0043231,GO:0043233,GO:0043280,GO:0043281,GO:0043299,GO:0043312,GO:0043412,GO:0043436,GO:0043900,GO:0043901,GO:0044093,GO:0044237,GO:0044238,GO:0044248,GO:0044260,GO:0044267,GO:0044281,GO:0044421,GO:0044422,GO:0044424,GO:0044433,GO:0044444,GO:0044446,GO:0044464,GO:0044703,GO:0045055,GO:0045137,GO:0045321,GO:0045454,GO:0045862,GO:0046425,GO:0046427,GO:0046546,GO:0046661,GO:0046903,GO:0046983,GO:0048232,GO:0048513,GO:0048518,GO:0048519,GO:0048522,GO:0048523,GO:0048583,GO:0048584,GO:0048608,GO:0048609,GO:0048731,GO:0048856,GO:0050789,GO:0050790,GO:0050794,GO:0050896,GO:0051171,GO:0051173,GO:0051179,GO:0051186,GO:0051187,GO:0051234,GO:0051239,GO:0051241,GO:0051246,GO:0051247,GO:0051336,GO:0051345,GO:0051604,GO:0051704,GO:0051707,GO:0051716,GO:0051920,GO:0052547,GO:0052548,GO:0055114,GO:0060205,GO:0060255,GO:0061458,GO:0065007,GO:0065008,GO:0065009,GO:0070013,GO:0070417,GO:0070887,GO:0071704,GO:0071840,GO:0072593,GO:0080090,GO:0097190,GO:0097237,GO:0097708,GO:0098542,GO:0098754,GO:0098869,GO:0099503,GO:0101002,GO:1901564,GO:1901605,GO:1902531,GO:1902533,GO:1904813,GO:1904892,GO:1904894,GO:1905936,GO:1905937,GO:1990748,GO:2000116,GO:2000241,GO:2000242,GO:2000254,GO:2000255,GO:2001056,GO:2001233,GO:2001235,GO:2001267,GO:2001269
## EC KEGG_ko
## 1 2.1.1.14 ko:K00549
## 2 1.11.1.15 ko:K03386
## KEGG_Pathway
## 1 ko00270,ko00450,ko01100,ko01110,ko01230,map00270,map00450,map01100,map01110,map01230
## 2 ko04214,map04214
## KEGG_Module KEGG_Reaction KEGG_rclass
## 1 M00017 R04405,R09365 RC00035,RC00113,RC01241
## 2 - - -
## BRITE KEGG_TC CAZy BiGG_Reaction
## 1 ko00000,ko00001,ko00002,ko01000 - - -
## 2 ko00000,ko00001,ko01000,ko04147 - - -
## PFAMs
## 1 Meth_synt_2
## 2 1-cysPrx_C,AhpC-TSAgogene <- merge(transcripts, eggnog, by=c("gene_id"), all=T)
gogene <- merge(gogene, kegg, by=c("gene_id"), all=T)
head(gogene,2)## gene_id seqnames
## 1 Pocillopora_acuta_HIv2___RNAseq.10002_t Pocillopora_acuta_HIv2___Sc0000013
## 2 Pocillopora_acuta_HIv2___RNAseq.10010_t Pocillopora_acuta_HIv2___Sc0000013
## start end width strand source type score.x phase
## 1 4542087 4551503 9417 + AUGUSTUS transcript NA NA
## 2 4639103 4647350 8248 + AUGUSTUS transcript NA NA
## ID
## 1 Pocillopora_acuta_HIv2___RNAseq.10002_t
## 2 Pocillopora_acuta_HIv2___RNAseq.10010_t
## transcript_id Parent seed_ortholog evalue
## 1 Pocillopora_acuta_HIv2___RNAseq.10002_t 45351.EDO27354 2.41e-93
## 2 Pocillopora_acuta_HIv2___RNAseq.10010_t 6087.XP_002166004.2 1.28e-38
## score.y
## 1 317
## 2 164
## eggNOG_OGs
## 1 COG0666@1|root,KOG0510@2759|Eukaryota,38G7Q@33154|Opisthokonta,3BCDU@33208|Metazoa
## 2 COG0666@1|root,KOG4177@2759|Eukaryota
## max_annot_lvl COG_category Description
## 1 33208|Metazoa DZ osmolarity-sensing cation channel activity
## 2 2759|Eukaryota I spectrin binding
## Preferred_name
## 1 TRPA1
## 2 -
## GOs
## 1 GO:0000302,GO:0001580,GO:0002791,GO:0002793,GO:0003008,GO:0003012,GO:0003674,GO:0004888,GO:0005034,GO:0005215,GO:0005216,GO:0005217,GO:0005244,GO:0005245,GO:0005261,GO:0005262,GO:0005488,GO:0005515,GO:0005575,GO:0005623,GO:0005886,GO:0005887,GO:0006810,GO:0006811,GO:0006812,GO:0006816,GO:0006873,GO:0006874,GO:0006875,GO:0006936,GO:0006939,GO:0006950,GO:0006979,GO:0007154,GO:0007165,GO:0007166,GO:0007204,GO:0007600,GO:0007602,GO:0007606,GO:0007610,GO:0007638,GO:0008150,GO:0008324,GO:0009266,GO:0009314,GO:0009408,GO:0009409,GO:0009410,GO:0009416,GO:0009453,GO:0009581,GO:0009582,GO:0009583,GO:0009593,GO:0009605,GO:0009612,GO:0009628,GO:0009636,GO:0009719,GO:0009966,GO:0009967,GO:0009987,GO:0010033,GO:0010035,GO:0010037,GO:0010243,GO:0010378,GO:0010646,GO:0010647,GO:0010817,GO:0014070,GO:0014074,GO:0014832,GO:0014848,GO:0015075,GO:0015085,GO:0015267,GO:0015276,GO:0015278,GO:0015318,GO:0016020,GO:0016021,GO:0016043,GO:0016048,GO:0016324,GO:0019233,GO:0019722,GO:0019725,GO:0019932,GO:0022607,GO:0022803,GO:0022832,GO:0022834,GO:0022836,GO:0022838,GO:0022839,GO:0022843,GO:0022857,GO:0022890,GO:0023041,GO:0023051,GO:0023052,GO:0023056,GO:0030001,GO:0030003,GO:0030424,GO:0031000,GO:0031224,GO:0031226,GO:0031644,GO:0031646,GO:0032024,GO:0032421,GO:0032501,GO:0032879,GO:0032880,GO:0032991,GO:0033554,GO:0033555,GO:0034220,GO:0034605,GO:0034702,GO:0034703,GO:0035556,GO:0035690,GO:0035774,GO:0036270,GO:0038023,GO:0040011,GO:0040040,GO:0042221,GO:0042330,GO:0042331,GO:0042391,GO:0042493,GO:0042542,GO:0042592,GO:0042752,GO:0042802,GO:0042995,GO:0043005,GO:0043052,GO:0043269,GO:0043270,GO:0043279,GO:0043933,GO:0044057,GO:0044070,GO:0044085,GO:0044425,GO:0044459,GO:0044464,GO:0045177,GO:0046677,GO:0046873,GO:0046883,GO:0046887,GO:0046957,GO:0048265,GO:0048518,GO:0048519,GO:0048522,GO:0048523,GO:0048583,GO:0048584,GO:0048878,GO:0050708,GO:0050714,GO:0050789,GO:0050794,GO:0050796,GO:0050801,GO:0050848,GO:0050850,GO:0050877,GO:0050896,GO:0050906,GO:0050907,GO:0050909,GO:0050912,GO:0050913,GO:0050951,GO:0050954,GO:0050955,GO:0050960,GO:0050961,GO:0050965,GO:0050966,GO:0050968,GO:0050974,GO:0050982,GO:0051046,GO:0051047,GO:0051049,GO:0051050,GO:0051179,GO:0051209,GO:0051222,GO:0051223,GO:0051234,GO:0051239,GO:0051240,GO:0051259,GO:0051260,GO:0051262,GO:0051282,GO:0051283,GO:0051289,GO:0051480,GO:0051606,GO:0051641,GO:0051649,GO:0051716,GO:0051930,GO:0051931,GO:0051969,GO:0052129,GO:0055065,GO:0055074,GO:0055080,GO:0055082,GO:0055085,GO:0060089,GO:0060341,GO:0060401,GO:0060402,GO:0061178,GO:0065003,GO:0065007,GO:0065008,GO:0070201,GO:0070417,GO:0070588,GO:0070838,GO:0070887,GO:0071241,GO:0071244,GO:0071310,GO:0071312,GO:0071313,GO:0071407,GO:0071415,GO:0071417,GO:0071466,GO:0071495,GO:0071840,GO:0071944,GO:0072347,GO:0072503,GO:0072507,GO:0072511,GO:0090087,GO:0090276,GO:0090277,GO:0097458,GO:0097553,GO:0097603,GO:0097604,GO:0098590,GO:0098655,GO:0098660,GO:0098662,GO:0098771,GO:0098796,GO:0098862,GO:0098900,GO:0098908,GO:0099094,GO:0099604,GO:0120025,GO:1901698,GO:1901699,GO:1901700,GO:1901701,GO:1902495,GO:1902531,GO:1902533,GO:1903522,GO:1903530,GO:1903532,GO:1903793,GO:1904058,GO:1904951,GO:1990351,GO:1990760
## 2 -
## EC KEGG_ko KEGG_Pathway KEGG_Module KEGG_Reaction KEGG_rclass
## 1 - ko:K04984 ko04750,map04750 - - -
## 2 - - - - - -
## BRITE KEGG_TC CAZy
## 1 ko00000,ko00001,ko04040 1.A.4.6.1,1.A.4.6.2,1.A.4.6.3,1.A.4.6.5 -
## 2 - - -
## BiGG_Reaction PFAMs KEGG_new
## 1 - Ank,Ank_2,Ank_3,Ank_4,Ion_trans K04984
## 2 - Ank,Ank_2,Ank_4,Ank_5,Ion_trans K04984geneInfo <- read.csv("../../output/WGCNA/WGCNA_ModuleMembership.csv") #this file was generated from the WGCNA analyses and contains the modules of interest
geneInfo<- plyr::rename(geneInfo, c("X"="gene_id"))
dim(geneInfo) # there are 9012 genes in our gene info file## [1] 9012 40geneInfo <- merge(gogene, geneInfo, by=c("gene_id"), all=T)Format GO terms to remove dashes and quotes and separate by semicolons (replace , with ;) in GOs column
geneInfo$GOs <- gsub(",", ";", geneInfo$GOs)
geneInfo$GOs <- gsub('"', "", geneInfo$GOs)
geneInfo$GOs <- gsub("-", NA, geneInfo$GOs)
geneInfo$KEGG_new[geneInfo$KEGG_new == ""] <- NAunique(geneInfo$moduleColor)## [1] "green" NA "blue" "salmon" "turquoise"
## [6] "yellow" "black" "red" "magenta" "lightcyan"
## [11] "purple" "brown" "pink" "midnightblue" "tan"
## [16] "cyan"geneInfo$Length<-lengths$transcript_lengths[match(geneInfo$gene_id, lengths$seqnames)]Frequency of GO terms
cal_biomin_terms <- read.csv("../../output/Biomineralization_goterms.csv")
head(cal_biomin_terms)## X.1 X GOterm over_represented_pvalue under_represented_pvalue numDEInCat
## 1 1 1 GO:0006325 0.009864096 1 2
## 2 2 2 GO:0016570 0.009864096 1 2
## 3 3 3 GO:0051276 0.009864096 1 2
## 4 4 4 GO:0060537 0.012664490 1 2
## 5 5 5 GO:0000278 0.012973357 1 2
## 6 6 6 GO:0007049 0.012973357 1 2
## numInCat term ontology bh_adjust
## 1 2 chromatin organization BP 0.8862339
## 2 2 histone modification BP 0.8862339
## 3 2 chromosome organization BP 0.8862339
## 4 2 muscle tissue development BP 0.8862339
## 5 2 mitotic cell cycle BP 0.8862339
## 6 2 cell cycle BP 0.8862339
## ParentTerm Factor
## 1 chromatin remodeling Biomin
## 2 macromolecule deacylation Biomin
## 3 chromatin remodeling Biomin
## 4 muscle tissue development Biomin
## 5 microtubule cytoskeleton organization involved in mitosis Biomin
## 6 microtubule cytoskeleton organization involved in mitosis Biomincal_up_terms <- read.csv("../../output/WGCNA/GO_analysis/goseq_pattern_calcification_filtered.csv")
cal_up_terms <- cal_up_terms %>%
mutate(Factor = "Up")
head(cal_up_terms)## X.1 X GOterm over_represented_pvalue under_represented_pvalue numDEInCat
## 1 1 1 GO:0000003 0 1 465
## 2 2 2 GO:0006139 0 1 608
## 3 3 3 GO:0006355 0 1 462
## 4 4 4 GO:0006725 0 1 650
## 5 5 5 GO:0006807 0 1 1215
## 6 6 6 GO:0006810 0 1 664
## numInCat term ontology bh_adjust
## 1 465 reproduction BP 0
## 2 608 nucleobase-containing compound metabolic process BP 0
## 3 462 regulation of DNA-templated transcription BP 0
## 4 650 cellular aromatic compound metabolic process BP 0
## 5 1215 nitrogen compound metabolic process BP 0
## 6 664 transport BP 0
## ParentTerm Factor
## 1 reproduction Up
## 2 metabolic process Up
## 3 regulation of biosynthetic process Up
## 4 metabolic process Up
## 5 metabolic process Up
## 6 localization Upcal_down_terms <-read.csv("../../output/WGCNA/GO_analysis/goseq_pattern_calcification_down_filtered.csv")
cal_down_terms<-cal_down_terms %>%
mutate(Factor = "Down")
head(cal_down_terms)## X.1 X GOterm over_represented_pvalue under_represented_pvalue numDEInCat
## 1 1 1 GO:0006139 0 1 485
## 2 2 2 GO:0006725 0 1 513
## 3 3 3 GO:0006807 0 1 910
## 4 4 4 GO:0006810 0 1 419
## 5 5 5 GO:0006950 0 1 364
## 6 6 6 GO:0006996 0 1 461
## numInCat term ontology bh_adjust
## 1 485 nucleobase-containing compound metabolic process BP 0
## 2 513 cellular aromatic compound metabolic process BP 0
## 3 910 nitrogen compound metabolic process BP 0
## 4 419 transport BP 0
## 5 364 response to stress BP 0
## 6 461 organelle organization BP 0
## ParentTerm Factor
## 1 gene expression Down
## 2 gene expression Down
## 3 gene expression Down
## 4 localization Down
## 5 response to stress Down
## 6 cellular component organization or biogenesis Downcolnames(cal_biomin_terms)## [1] "X.1" "X"
## [3] "GOterm" "over_represented_pvalue"
## [5] "under_represented_pvalue" "numDEInCat"
## [7] "numInCat" "term"
## [9] "ontology" "bh_adjust"
## [11] "ParentTerm" "Factor"colnames(cal_up_terms)## [1] "X.1" "X"
## [3] "GOterm" "over_represented_pvalue"
## [5] "under_represented_pvalue" "numDEInCat"
## [7] "numInCat" "term"
## [9] "ontology" "bh_adjust"
## [11] "ParentTerm" "Factor"colnames(cal_down_terms)## [1] "X.1" "X"
## [3] "GOterm" "over_represented_pvalue"
## [5] "under_represented_pvalue" "numDEInCat"
## [7] "numInCat" "term"
## [9] "ontology" "bh_adjust"
## [11] "ParentTerm" "Factor"Merge biomineralization, up and down-regulation of calcification GOterms
all_terms<- merge(cal_up_terms,cal_down_terms, by=c("Factor","GOterm","X.1","X","GOterm","over_represented_pvalue","under_represented_pvalue","numDEInCat","numInCat","term","ontology","bh_adjust","ParentTerm"),all=T)
all_terms<- merge(all_terms,cal_biomin_terms, by=c("Factor","GOterm","X.1","X","GOterm","over_represented_pvalue","under_represented_pvalue","numDEInCat","numInCat","term","ontology","bh_adjust","ParentTerm"),all=T)
all_terms$GOterm<-as.factor(all_terms$GOterm)
head(all_terms)## Factor GOterm X.1 X over_represented_pvalue under_represented_pvalue
## 1 Biomin GO:0000003 300 343 0.54128075 0.8456150
## 2 Biomin GO:0000041 417 565 1.00000000 0.8535542
## 3 Biomin GO:0000122 65 70 0.12290664 1.0000000
## 4 Biomin GO:0000132 107 122 0.15108190 1.0000000
## 5 Biomin GO:0000226 256 298 0.37841205 0.9418666
## 6 Biomin GO:0000278 5 5 0.01297336 1.0000000
## numDEInCat numInCat term
## 1 1 5 reproduction
## 2 0 1 transition metal ion transport
## 3 1 1 negative regulation of transcription by RNA polymerase II
## 4 1 1 establishment of mitotic spindle orientation
## 5 1 3 microtubule cytoskeleton organization
## 6 2 2 mitotic cell cycle
## ontology bh_adjust ParentTerm
## 1 BP 1.0000000 female sex differentiation
## 2 BP 1.0000000 calcium ion import
## 3 BP 0.8862339 negative regulation of biological process
## 4 BP 0.8862339 establishment of mitotic spindle orientation
## 5 BP 1.0000000 microtubule-based process
## 6 BP 0.8862339 microtubule cytoskeleton organization involved in mitosistail(all_terms)## Factor GOterm X.1 X over_represented_pvalue
## 4614 Up GO:2001242 922 983 5.023572e-31
## 4615 Up GO:2001243 1342 1500 2.171774e-18
## 4616 Up GO:2001251 978 1045 3.318296e-29
## 4617 Up GO:2001252 713 759 1.527374e-42
## 4618 Up GO:2001257 1514 1737 2.346995e-15
## 4619 Up GO:2001259 2050 2535 2.136874e-09
## under_represented_pvalue numDEInCat numInCat
## 4614 1 42 42
## 4615 1 24 24
## 4616 1 43 43
## 4617 1 61 61
## 4618 1 21 21
## 4619 1 12 12
## term ontology
## 4614 regulation of intrinsic apoptotic signaling pathway BP
## 4615 negative regulation of intrinsic apoptotic signaling pathway BP
## 4616 negative regulation of chromosome organization BP
## 4617 positive regulation of chromosome organization BP
## 4618 regulation of cation channel activity BP
## 4619 positive regulation of cation channel activity BP
## bh_adjust ParentTerm
## 4614 5.462041e-30 regulation of cell death
## 4615 1.558880e-17 regulation of cell death
## 4616 3.389807e-28 regulation of cellular component organization
## 4617 2.152863e-41 regulation of cellular component organization
## 4618 1.467562e-14 regulation of localization
## 4619 9.175464e-09 regulation of localizationgoterms_shared <- all_terms %>%
group_by(GOterm) %>%
dplyr::summarise(
ParentTerm = paste(unique(ParentTerm), collapse = ", "),
Factor = paste(unique(Factor), collapse = ", ")
)
dim(goterms_shared)## [1] 2322 3write.csv(goterms_shared, "../../output/WGCNA/GO_analysis/Merged_GOterms_factor_ParentTerm.csv")# How is this file made?
cal_freq_terms <-read.csv("../../output/WGCNA/Kristen_Old/goseq_pattern_calcification_all.csv")
head(cal_freq_terms)## ParentTerm Calcification.direction
## 1 negative regulation of organelle organization Up
## 2 negative regulation of protein modification process Up
## 3 anion transport Up
## 4 sensory organ morphogenesis Up
## 5 cytokine-mediated signaling pathway Up
## 6 biological regulation Up
## Number.of.terms Frequency.of.shared.biomin.GOterms
## 1 39 12
## 2 26 11
## 3 23 6
## 4 21 3
## 5 20 2
## 6 19 9
## Proportion.of.shared.GO.terms.with.biomineralization.genes
## 1 0.3076923
## 2 0.4230769
## 3 0.2608696
## 4 0.1428571
## 5 0.1000000
## 6 0.4736842
## Percentage.of.shared.GO.terms.with.biomineralization.genes
## 1 30.76923
## 2 42.30769
## 3 26.08696
## 4 14.28571
## 5 10.00000
## 6 47.36842##Frequency >10 GO terms upregulation
cal_freq_terms_filtered_up <- cal_freq_terms %>%
filter(Number.of.terms>=10) %>%
filter(Calcification.direction=="Up")
cal_freq_terms_filtered_up## ParentTerm
## 1 negative regulation of organelle organization
## 2 negative regulation of protein modification process
## 3 anion transport
## 4 sensory organ morphogenesis
## 5 cytokine-mediated signaling pathway
## 6 biological regulation
## 7 establishment of protein localization to organelle
## 8 positive regulation of phosphate metabolic process
## 9 proteasomal protein catabolic process
## 10 cellular metal ion homeostasis
## 11 negative regulation of transcription by RNA polymerase II
## 12 post-embryonic animal morphogenesis
## 13 regulation of intracellular protein transport
## 14 regulation of protein kinase activity
## 15 negative regulation of intracellular signal transduction
## 16 sensory perception
## 17 transcription by RNA polymerase II
## 18 programmed cell death
## 19 meiotic cell cycle process
## 20 muscle structure development
## 21 organic acid biosynthetic process
## 22 ncRNA metabolic process
## Calcification.direction Number.of.terms Frequency.of.shared.biomin.GOterms
## 1 Up 39 12
## 2 Up 26 11
## 3 Up 23 6
## 4 Up 21 3
## 5 Up 20 2
## 6 Up 19 9
## 7 Up 18 5
## 8 Up 18 5
## 9 Up 18 3
## 10 Up 17 17
## 11 Up 17 17
## 12 Up 17 17
## 13 Up 17 5
## 14 Up 16 15
## 15 Up 15 15
## 16 Up 14 14
## 17 Up 14 1
## 18 Up 12 5
## 19 Up 11 6
## 20 Up 11 10
## 21 Up 11 7
## 22 Up 10 3
## Proportion.of.shared.GO.terms.with.biomineralization.genes
## 1 0.30769231
## 2 0.42307692
## 3 0.26086957
## 4 0.14285714
## 5 0.10000000
## 6 0.47368421
## 7 0.27777778
## 8 0.27777778
## 9 0.16666667
## 10 1.00000000
## 11 1.00000000
## 12 1.00000000
## 13 0.29411765
## 14 0.93750000
## 15 1.00000000
## 16 1.00000000
## 17 0.07142857
## 18 0.41666667
## 19 0.54545455
## 20 0.90909091
## 21 0.63636364
## 22 0.30000000
## Percentage.of.shared.GO.terms.with.biomineralization.genes
## 1 30.769231
## 2 42.307692
## 3 26.086957
## 4 14.285714
## 5 10.000000
## 6 47.368421
## 7 27.777778
## 8 27.777778
## 9 16.666667
## 10 100.000000
## 11 100.000000
## 12 100.000000
## 13 29.411765
## 14 93.750000
## 15 100.000000
## 16 100.000000
## 17 7.142857
## 18 41.666667
## 19 54.545455
## 20 90.909091
## 21 63.636364
## 22 30.000000Figure
#counts$Direction.of.flat.origin<- factor(counts$Direction.of.flat.origin, levels =c("up","no pattern","down"))
#counts$Module<- factor(counts$Module, levels=c("Blue","Brown","Greenyellow","Cyan","Pink","Magenta","Lightcyan","Midnight blue","Purple","Turquiose","Red","Black"))
freq_fig_up<-ggplot(cal_freq_terms_filtered_up, aes(y=Number.of.terms,x=reorder(ParentTerm, Number.of.terms), fill=Percentage.of.shared.GO.terms.with.biomineralization.genes,group=1))+
#facet_wrap(~Calcification.direction, nrow = 1)+
geom_point(size=5, alpha=1, pch=21,color="black")+
geom_segment(aes(x=ParentTerm, xend=ParentTerm, y=0, yend=Number.of.terms)) +
geom_hline(yintercept = 0, linetype="solid", color = 'black', size=0.5, show.legend = TRUE)+
coord_flip()+
scale_y_continuous(expression(GO~term~counts),limits=c(0,40))+
#scale_color_manual("Direction.of.flat.origin",values= c("up"="#b2182b","no pattern"="grey","down" ="#67a9cf"))+
#scale_fill_manual("Direction.of.flat.origin",values= c("up"="#b2182b","no pattern"="grey","down" ="#67a9cf"))+
scale_fill_gradientn(colours=c("white","#fddbc7","#f4a582","#d6604d","#b2182b"), na.value = "grey98",limits = c(0, 100))+
#scale_color_gradientn(colours=c("#b2182b","#fddbc7","white","#d1e5f0","#67a9cf", "#67a9cf", "#2166ac"), na.value = "grey98",limits = c(-0, 40))+
theme_classic()+
theme(axis.text.x=element_text(vjust=0.5, hjust=0.95,size=12),
plot.title = element_text(margin = margin(t = 10, b = 10), hjust=0.5),
panel.background= element_rect(fill=NA, color='black'),
legend.title = element_blank(),
axis.text.y = element_text(vjust=0.5, size=12), #making the axis text larger
axis.title.x = element_blank(),#making the axis title larger
axis.title.y = element_blank(),
strip.text = element_text(size=12))#making the axis title larger ## Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
## ℹ Please use `linewidth` instead.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.freq_fig_up
ggsave(filename="../../output/WGCNA/GO_analysis/freq_fig_up_extra_filtered.png", plot=freq_fig_up, dpi=300, height=6, units="in", limitsize=FALSE)## Saving 7 x 6 in imageFrequency >10 GO terms downregulation
cal_freq_terms_filtered_down <- cal_freq_terms %>%
filter(Number.of.terms>=10) %>%
filter(Calcification.direction=="Down")
cal_freq_terms_filtered_down## ParentTerm
## 1 aging
## 2 biological regulation
## 3 cation transport
## 4 development of primary sexual characteristics
## 5 embryonic organ development
## 6 histone modification
## 7 negative regulation of catabolic process
## 8 negative regulation of cell development
## 9 negative regulation of mitotic cell cycle
## 10 negative regulation of neurogenesis
## 11 negative regulation of transcription by RNA polymerase II
## 12 positive regulation of kinase activity
## 13 positive regulation of multi-organism process
## 14 post-embryonic animal morphogenesis
## 15 protein localization to organelle
## 16 purine ribonucleotide metabolic process
## 17 regulation of cell size
## 18 regulation of gene expression, epigenetic
## 19 regulation of ion transport
## 20 regulation of neuron death
## 21 response to xenobiotic stimulus
## 22 ribose phosphate metabolic process
## 23 transcription by RNA polymerase II
## 24 transmembrane receptor protein tyrosine kinase signaling pathway
## Calcification.direction Number.of.terms Frequency.of.shared.biomin.GOterms
## 1 Down 12 11
## 2 Down 16 9
## 3 Down 14 14
## 4 Down 10 8
## 5 Down 19 19
## 6 Down 39 27
## 7 Down 13 12
## 8 Down 10 10
## 9 Down 10 7
## 10 Down 13 11
## 11 Down 17 17
## 12 Down 20 17
## 13 Down 10 6
## 14 Down 19 19
## 15 Down 14 9
## 16 Down 13 7
## 17 Down 20 18
## 18 Down 27 26
## 19 Down 17 7
## 20 Down 13 5
## 21 Down 25 19
## 22 Down 22 12
## 23 Down 14 1
## 24 Down 23 21
## Proportion.of.shared.GO.terms.with.biomineralization.genes
## 1 0.91666667
## 2 0.56250000
## 3 1.00000000
## 4 0.80000000
## 5 1.00000000
## 6 0.69230769
## 7 0.92307692
## 8 1.00000000
## 9 0.70000000
## 10 0.84615385
## 11 1.00000000
## 12 0.85000000
## 13 0.60000000
## 14 1.00000000
## 15 0.64285714
## 16 0.53846154
## 17 0.90000000
## 18 0.96296296
## 19 0.41176471
## 20 0.38461538
## 21 0.76000000
## 22 0.54545455
## 23 0.07142857
## 24 0.91304348
## Percentage.of.shared.GO.terms.with.biomineralization.genes
## 1 91.666667
## 2 56.250000
## 3 100.000000
## 4 80.000000
## 5 100.000000
## 6 69.230769
## 7 92.307692
## 8 100.000000
## 9 70.000000
## 10 84.615385
## 11 100.000000
## 12 85.000000
## 13 60.000000
## 14 100.000000
## 15 64.285714
## 16 53.846154
## 17 90.000000
## 18 96.296296
## 19 41.176471
## 20 38.461538
## 21 76.000000
## 22 54.545455
## 23 7.142857
## 24 91.304348Figure
freq_fig_down<-ggplot(cal_freq_terms_filtered_down, aes(y=Number.of.terms,x=reorder(ParentTerm, Number.of.terms), fill=Percentage.of.shared.GO.terms.with.biomineralization.genes))+
#facet_wrap(~Calcification.direction, nrow = 1)+
geom_point(size=5, alpha=1, pch=21,color="black")+
geom_segment(aes(x=ParentTerm, xend=ParentTerm, y=0, yend=Number.of.terms)) +
#geom_hline(yintercept = 0, linetype="solid", color = 'black', size=0.5, show.legend = TRUE)+
coord_flip()+
scale_y_continuous(expression(GO~term~counts),limits=c(0,40))+
#scale_color_manual("Direction.of.flat.origin",values= c("up"="#b2182b","no pattern"="grey","down" ="#67a9cf"))+
#scale_fill_manual("Direction.of.flat.origin",values= c("up"="#b2182b","no pattern"="grey","down" ="#67a9cf"))+
scale_fill_gradientn(colours=c("white","#d1e5f0","#92c5de","#4393c3","#2166ac"), na.value = "grey98",limits = c(0, 100))+
#scale_fill_gradientn(colours=c("#fddbc7","#f4a582","#d6604d","#b2182b"), na.value = "grey98",limits = c(10, 40))+
#scale_color_gradientn(colours=c("#b2182b","#fddbc7","white","#d1e5f0","#67a9cf", "#67a9cf", "#2166ac"), na.value = "grey98",limits = c(-0, 40))+
theme_classic()+
theme(axis.text.x=element_text(vjust=0.5, hjust=0.95, size=12),
plot.title = element_text(margin = margin(t = 10, b = 10), hjust=0.5),
panel.background= element_rect(fill=NA, color='black'),
legend.title = element_blank(),
axis.text.y = element_text(vjust=0.5, size=12), #making the axis text larger
axis.title.x = element_blank(),#making the axis title larger
axis.title.y = element_blank(),
strip.text = element_text(size=12))#making the axis title larger
freq_fig_down
ggsave(filename="../../output/WGCNA/GO_analysis/freq_fig_down_extra_filtered.png", plot=freq_fig_down, dpi=300, height=6, units="in", limitsize=FALSE)## Saving 7 x 6 in imagecompare_figs<-cowplot::plot_grid(freq_fig_up, freq_fig_down, nrow=2, align="v")
compare_figs
##Frequency all GO terms upregulation
cal_freq_terms_filtered_up_all <- cal_freq_terms %>%
filter(Calcification.direction=="Up")
cal_freq_terms_filtered_up_all## ParentTerm
## 1 negative regulation of organelle organization
## 2 negative regulation of protein modification process
## 3 anion transport
## 4 sensory organ morphogenesis
## 5 cytokine-mediated signaling pathway
## 6 biological regulation
## 7 establishment of protein localization to organelle
## 8 positive regulation of phosphate metabolic process
## 9 proteasomal protein catabolic process
## 10 cellular metal ion homeostasis
## 11 negative regulation of transcription by RNA polymerase II
## 12 post-embryonic animal morphogenesis
## 13 regulation of intracellular protein transport
## 14 regulation of protein kinase activity
## 15 negative regulation of intracellular signal transduction
## 16 sensory perception
## 17 transcription by RNA polymerase II
## 18 programmed cell death
## 19 meiotic cell cycle process
## 20 muscle structure development
## 21 organic acid biosynthetic process
## 22 ncRNA metabolic process
## 23 positive regulation of cell development
## 24 oogenesis
## 25 sex differentiation
## 26 cell-cell signaling
## 27 defense response
## 28 protein ubiquitination
## 29 regulation of cell migration
## 30 cell fate commitment
## 31 epithelial cell development
## 32 metabolic process
## 33 nitrogen compound metabolic process
## 34 regulation of cellular response to stress
## 35 regulation of immune response
## 36 chemotaxis
## 37 positive regulation of growth
## 38 purine ribonucleotide metabolic process
## 39 response to other organism
## 40 actin cytoskeleton organization
## 41 lipid biosynthetic process
## 42 response to extracellular stimulus
## 43 response to radiation
## 44 amide biosynthetic process
## 45 cell population proliferation
## 46 localization of cell
## 47 response to endogenous stimulus
## 48 tube development
## 49 biological adhesion
## 50 brain development
## 51 carbohydrate derivative biosynthetic process
## 52 cell activation
## 53 gene expression
## 54 generation of precursor metabolites and energy
## 55 microtubule-based process
## 56 proteolysis
## 57 behavior
## 58 carbohydrate metabolic process
## 59 cell division
## 60 cellular process
## 61 developmental maturation
## 62 drug metabolic process
## 63 import into cell
## 64 organic hydroxy compound metabolic process
## 65 oxidation-reduction process
## 66 sulfur compound metabolic process
## Calcification.direction Number.of.terms Frequency.of.shared.biomin.GOterms
## 1 Up 39 12
## 2 Up 26 11
## 3 Up 23 6
## 4 Up 21 3
## 5 Up 20 2
## 6 Up 19 9
## 7 Up 18 5
## 8 Up 18 5
## 9 Up 18 3
## 10 Up 17 17
## 11 Up 17 17
## 12 Up 17 17
## 13 Up 17 5
## 14 Up 16 15
## 15 Up 15 15
## 16 Up 14 14
## 17 Up 14 1
## 18 Up 12 5
## 19 Up 11 6
## 20 Up 11 10
## 21 Up 11 7
## 22 Up 10 3
## 23 Up 9 9
## 24 Up 8 6
## 25 Up 8 6
## 26 Up 7 7
## 27 Up 7 4
## 28 Up 7 4
## 29 Up 7 7
## 30 Up 6 4
## 31 Up 6 6
## 32 Up 6 6
## 33 Up 6 5
## 34 Up 6 4
## 35 Up 6 1
## 36 Up 5 5
## 37 Up 5 5
## 38 Up 5 5
## 39 Up 5 2
## 40 Up 4 2
## 41 Up 4 2
## 42 Up 4 4
## 43 Up 4 3
## 44 Up 3 2
## 45 Up 3 1
## 46 Up 3 3
## 47 Up 3 3
## 48 Up 3 3
## 49 Up 2 2
## 50 Up 2 2
## 51 Up 2 1
## 52 Up 2 0
## 53 Up 2 1
## 54 Up 2 0
## 55 Up 2 2
## 56 Up 2 2
## 57 Up 1 1
## 58 Up 1 0
## 59 Up 1 1
## 60 Up 1 1
## 61 Up 1 1
## 62 Up 1 1
## 63 Up 1 1
## 64 Up 1 1
## 65 Up 1 1
## 66 Up 1 1
## Proportion.of.shared.GO.terms.with.biomineralization.genes
## 1 0.30769231
## 2 0.42307692
## 3 0.26086957
## 4 0.14285714
## 5 0.10000000
## 6 0.47368421
## 7 0.27777778
## 8 0.27777778
## 9 0.16666667
## 10 1.00000000
## 11 1.00000000
## 12 1.00000000
## 13 0.29411765
## 14 0.93750000
## 15 1.00000000
## 16 1.00000000
## 17 0.07142857
## 18 0.41666667
## 19 0.54545455
## 20 0.90909091
## 21 0.63636364
## 22 0.30000000
## 23 1.00000000
## 24 0.75000000
## 25 0.75000000
## 26 1.00000000
## 27 0.57142857
## 28 0.57142857
## 29 1.00000000
## 30 0.66666667
## 31 1.00000000
## 32 1.00000000
## 33 0.83333333
## 34 0.66666667
## 35 0.16666667
## 36 1.00000000
## 37 1.00000000
## 38 1.00000000
## 39 0.40000000
## 40 0.50000000
## 41 0.50000000
## 42 1.00000000
## 43 0.75000000
## 44 0.66666667
## 45 0.33333333
## 46 1.00000000
## 47 1.00000000
## 48 1.00000000
## 49 1.00000000
## 50 1.00000000
## 51 0.50000000
## 52 0.00000000
## 53 0.50000000
## 54 0.00000000
## 55 1.00000000
## 56 0.50000000
## 57 1.00000000
## 58 0.00000000
## 59 1.00000000
## 60 1.00000000
## 61 1.00000000
## 62 1.00000000
## 63 1.00000000
## 64 1.00000000
## 65 1.00000000
## 66 1.00000000
## Percentage.of.shared.GO.terms.with.biomineralization.genes
## 1 30.769231
## 2 42.307692
## 3 26.086957
## 4 14.285714
## 5 10.000000
## 6 47.368421
## 7 27.777778
## 8 27.777778
## 9 16.666667
## 10 100.000000
## 11 100.000000
## 12 100.000000
## 13 29.411765
## 14 93.750000
## 15 100.000000
## 16 100.000000
## 17 7.142857
## 18 41.666667
## 19 54.545455
## 20 90.909091
## 21 63.636364
## 22 30.000000
## 23 100.000000
## 24 75.000000
## 25 75.000000
## 26 100.000000
## 27 57.142857
## 28 57.142857
## 29 100.000000
## 30 66.666667
## 31 100.000000
## 32 100.000000
## 33 83.333333
## 34 66.666667
## 35 16.666667
## 36 100.000000
## 37 100.000000
## 38 100.000000
## 39 40.000000
## 40 50.000000
## 41 50.000000
## 42 100.000000
## 43 75.000000
## 44 66.666667
## 45 33.333333
## 46 100.000000
## 47 100.000000
## 48 100.000000
## 49 100.000000
## 50 100.000000
## 51 50.000000
## 52 0.000000
## 53 50.000000
## 54 0.000000
## 55 100.000000
## 56 50.000000
## 57 100.000000
## 58 0.000000
## 59 100.000000
## 60 100.000000
## 61 100.000000
## 62 100.000000
## 63 100.000000
## 64 100.000000
## 65 100.000000
## 66 100.000000###Figure
#counts$Direction.of.flat.origin<- factor(counts$Direction.of.flat.origin, levels =c("up","no pattern","down"))
#counts$Module<- factor(counts$Module, levels=c("Blue","Brown","Greenyellow","Cyan","Pink","Magenta","Lightcyan","Midnight blue","Purple","Turquiose","Red","Black"))
freq_fig_up<-ggplot(cal_freq_terms_filtered_up_all, aes(y=Number.of.terms,x=reorder(ParentTerm, Number.of.terms), fill=Percentage.of.shared.GO.terms.with.biomineralization.genes,group=1))+
#facet_wrap(~Calcification.direction, nrow = 1)+
geom_point(size=5, alpha=1, pch=21,color="black")+
geom_segment(aes(x=ParentTerm, xend=ParentTerm, y=0, yend=Number.of.terms)) +
geom_hline(yintercept = 0, linetype="solid", color = 'black', size=0.5, show.legend = TRUE)+
coord_flip()+
scale_y_continuous(expression(GO~term~counts),limits=c(0,40))+
#scale_color_manual("Direction.of.flat.origin",values= c("up"="#b2182b","no pattern"="grey","down" ="#67a9cf"))+
#scale_fill_manual("Direction.of.flat.origin",values= c("up"="#b2182b","no pattern"="grey","down" ="#67a9cf"))+
scale_fill_gradientn(colours=c("white","#fddbc7","#f4a582","#d6604d","#b2182b"), na.value = "grey98",limits = c(0, 100))+
#scale_color_gradientn(colours=c("#b2182b","#fddbc7","white","#d1e5f0","#67a9cf", "#67a9cf", "#2166ac"), na.value = "grey98",limits = c(-0, 40))+
theme_classic()+
theme(axis.text.x=element_text(vjust=0.5, hjust=0.95,size=12),
plot.title = element_text(margin = margin(t = 10, b = 10), hjust=0.5),
panel.background= element_rect(fill=NA, color='black'),
legend.title = element_blank(),
axis.text.y = element_text(vjust=0.5, size=12), #making the axis text larger
axis.title.x = element_blank(),#making the axis title larger
axis.title.y = element_blank(),
strip.text = element_text(size=12))#making the axis title larger
freq_fig_up
##Frequency all GO terms downregulation
cal_freq_terms_filtered_down_all <- cal_freq_terms %>%
filter(Calcification.direction=="Down")
cal_freq_terms_filtered_down_all## ParentTerm
## 1 aging
## 2 axon guidance
## 3 biological regulation
## 4 carbohydrate derivative biosynthetic process
## 5 carbohydrate metabolic process
## 6 cation transport
## 7 cellular lipid metabolic process
## 8 cellular process
## 9 cellular response to abiotic stimulus
## 10 cellular response to external stimulus
## 11 cellular response to hormone stimulus
## 12 detection of stimulus
## 13 development of primary sexual characteristics
## 14 developmental growth
## 15 DNA repair
## 16 drug metabolic process
## 17 embryonic organ development
## 18 establishment or maintenance of cell polarity
## 19 exocytosis
## 20 gastrulation
## 21 head development
## 22 histone modification
## 23 import into cell
## 24 leukocyte activation
## 25 locomotory behavior
## 26 macromolecule metabolic process
## 27 memory
## 28 metabolic process
## 29 morphogenesis of embryonic epithelium
## 30 negative regulation of catabolic process
## 31 negative regulation of cell development
## 32 negative regulation of mitotic cell cycle
## 33 negative regulation of neurogenesis
## 34 negative regulation of transcription by RNA polymerase II
## 35 nucleotide metabolic process
## 36 organonitrogen compound metabolic process
## 37 oxidation-reduction process
## 38 positive regulation of cell motility
## 39 positive regulation of defense response
## 40 positive regulation of immune response
## 41 positive regulation of kinase activity
## 42 positive regulation of multi-organism process
## 43 post-embryonic animal morphogenesis
## 44 protein localization to organelle
## 45 protein modification by small protein conjugation or removal
## 46 purine ribonucleotide metabolic process
## 47 regulation of actin cytoskeleton organization
## 48 regulation of cell population proliferation
## 49 regulation of cell size
## 50 regulation of cell-cell adhesion
## 51 regulation of gene expression, epigenetic
## 52 regulation of ion transport
## 53 regulation of microtubule-based process
## 54 regulation of neuron death
## 55 regulation of peptide secretion
## 56 regulation of Ras protein signal transduction
## 57 response to ionizing radiation
## 58 response to xenobiotic stimulus
## 59 ribose phosphate metabolic process
## 60 rRNA metabolic process
## 61 transcription by RNA polymerase II
## 62 transmembrane receptor protein tyrosine kinase signaling pathway
## 63 tube development
## Calcification.direction Number.of.terms Frequency.of.shared.biomin.GOterms
## 1 Down 12 11
## 2 Down 7 7
## 3 Down 16 9
## 4 Down 2 1
## 5 Down 1 0
## 6 Down 14 14
## 7 Down 3 1
## 8 Down 1 1
## 9 Down 2 0
## 10 Down 5 5
## 11 Down 4 4
## 12 Down 6 3
## 13 Down 10 8
## 14 Down 8 6
## 15 Down 8 5
## 16 Down 1 1
## 17 Down 19 19
## 18 Down 1 1
## 19 Down 9 1
## 20 Down 7 7
## 21 Down 2 2
## 22 Down 39 27
## 23 Down 1 1
## 24 Down 2 0
## 25 Down 2 2
## 26 Down 3 3
## 27 Down 9 5
## 28 Down 6 6
## 29 Down 7 7
## 30 Down 13 12
## 31 Down 10 10
## 32 Down 10 7
## 33 Down 13 11
## 34 Down 17 17
## 35 Down 7 7
## 36 Down 6 4
## 37 Down 1 1
## 38 Down 8 7
## 39 Down 7 4
## 40 Down 9 1
## 41 Down 20 17
## 42 Down 10 6
## 43 Down 19 19
## 44 Down 14 9
## 45 Down 6 4
## 46 Down 13 7
## 47 Down 4 2
## 48 Down 4 1
## 49 Down 20 18
## 50 Down 2 0
## 51 Down 27 26
## 52 Down 17 7
## 53 Down 3 2
## 54 Down 13 5
## 55 Down 8 0
## 56 Down 7 6
## 57 Down 7 5
## 58 Down 25 19
## 59 Down 22 12
## 60 Down 6 1
## 61 Down 14 1
## 62 Down 23 21
## 63 Down 3 3
## Proportion.of.shared.GO.terms.with.biomineralization.genes
## 1 0.91666667
## 2 1.00000000
## 3 0.56250000
## 4 0.50000000
## 5 0.00000000
## 6 1.00000000
## 7 0.33333333
## 8 1.00000000
## 9 0.00000000
## 10 1.00000000
## 11 1.00000000
## 12 0.50000000
## 13 0.80000000
## 14 0.75000000
## 15 0.62500000
## 16 1.00000000
## 17 1.00000000
## 18 1.00000000
## 19 0.11111111
## 20 1.00000000
## 21 1.00000000
## 22 0.69230769
## 23 1.00000000
## 24 0.00000000
## 25 1.00000000
## 26 1.00000000
## 27 0.55555556
## 28 1.00000000
## 29 1.00000000
## 30 0.92307692
## 31 1.00000000
## 32 0.70000000
## 33 0.84615385
## 34 1.00000000
## 35 1.00000000
## 36 0.66666667
## 37 1.00000000
## 38 0.87500000
## 39 0.57142857
## 40 0.11111111
## 41 0.85000000
## 42 0.60000000
## 43 1.00000000
## 44 0.64285714
## 45 0.66666667
## 46 0.53846154
## 47 0.50000000
## 48 0.25000000
## 49 0.90000000
## 50 0.00000000
## 51 0.96296296
## 52 0.41176471
## 53 0.66666667
## 54 0.38461538
## 55 0.00000000
## 56 0.85714286
## 57 0.71428571
## 58 0.76000000
## 59 0.54545455
## 60 0.16666667
## 61 0.07142857
## 62 0.91304348
## 63 1.00000000
## Percentage.of.shared.GO.terms.with.biomineralization.genes
## 1 91.666667
## 2 100.000000
## 3 56.250000
## 4 50.000000
## 5 0.000000
## 6 100.000000
## 7 33.333333
## 8 100.000000
## 9 0.000000
## 10 100.000000
## 11 100.000000
## 12 50.000000
## 13 80.000000
## 14 75.000000
## 15 62.500000
## 16 100.000000
## 17 100.000000
## 18 100.000000
## 19 11.111111
## 20 100.000000
## 21 100.000000
## 22 69.230769
## 23 100.000000
## 24 0.000000
## 25 100.000000
## 26 100.000000
## 27 55.555556
## 28 100.000000
## 29 100.000000
## 30 92.307692
## 31 100.000000
## 32 70.000000
## 33 84.615385
## 34 100.000000
## 35 100.000000
## 36 66.666667
## 37 100.000000
## 38 87.500000
## 39 57.142857
## 40 11.111111
## 41 85.000000
## 42 60.000000
## 43 100.000000
## 44 64.285714
## 45 66.666667
## 46 53.846154
## 47 50.000000
## 48 25.000000
## 49 90.000000
## 50 0.000000
## 51 96.296296
## 52 41.176471
## 53 66.666667
## 54 38.461538
## 55 0.000000
## 56 85.714286
## 57 71.428571
## 58 76.000000
## 59 54.545455
## 60 16.666667
## 61 7.142857
## 62 91.304348
## 63 100.000000###Figure
freq_fig_down<-ggplot(cal_freq_terms_filtered_down_all, aes(y=Number.of.terms,x=reorder(ParentTerm, Number.of.terms), fill=Percentage.of.shared.GO.terms.with.biomineralization.genes))+
#facet_wrap(~Calcification.direction, nrow = 1)+
geom_point(size=5, alpha=1, pch=21,color="black")+
geom_segment(aes(x=ParentTerm, xend=ParentTerm, y=0, yend=Number.of.terms)) +
#geom_hline(yintercept = 0, linetype="solid", color = 'black', size=0.5, show.legend = TRUE)+
coord_flip()+
scale_y_continuous(expression(GO~term~counts),limits=c(0,40))+
#scale_color_manual("Direction.of.flat.origin",values= c("up"="#b2182b","no pattern"="grey","down" ="#67a9cf"))+
#scale_fill_manual("Direction.of.flat.origin",values= c("up"="#b2182b","no pattern"="grey","down" ="#67a9cf"))+
scale_fill_gradientn(colours=c("white","#d1e5f0","#92c5de","#4393c3","#2166ac"), na.value = "grey98",limits = c(0, 100))+
#scale_fill_gradientn(colours=c("#fddbc7","#f4a582","#d6604d","#b2182b"), na.value = "grey98",limits = c(10, 40))+
#scale_color_gradientn(colours=c("#b2182b","#fddbc7","white","#d1e5f0","#67a9cf", "#67a9cf", "#2166ac"), na.value = "grey98",limits = c(-0, 40))+
theme_classic()+
theme(axis.text.x=element_text(vjust=0.5, hjust=0.95, size=12),
plot.title = element_text(margin = margin(t = 10, b = 10), hjust=0.5),
panel.background= element_rect(fill=NA, color='black'),
legend.title = element_blank(),
axis.text.y = element_text(vjust=0.5, size=12), #making the axis text larger
axis.title.x = element_blank(),#making the axis title larger
axis.title.y = element_blank(),
strip.text = element_text(size=12))#making the axis title larger
freq_fig_down
compare_figs_all<-cowplot::plot_grid(freq_fig_up, freq_fig_down, ncol=2, align="h")
compare_figs_all
Biomineralization toolkit present in modules
biomin <-read.csv("../../output/Biomin_blast_Pocillopora_acuta_best_hit.csv")wgcnamod <-read.csv("../../output/WGCNA/WGCNA_ModuleMembership.csv")
wgcnamod <- plyr::rename(wgcnamod, c("X"="Pocillopora_acuta_best_hit"))biomin_mod <- merge(biomin, wgcnamod, by=c("Pocillopora_acuta_best_hit"), all=F)
head(biomin_mod)## Pocillopora_acuta_best_hit accessionnumber.geneID
## 1 Pocillopora_acuta_HIv2___RNAseq.g10093.t2 XP_022804785.1
## 2 Pocillopora_acuta_HIv2___RNAseq.g11609.t1 P33_g8985
## 3 Pocillopora_acuta_HIv2___RNAseq.g13172.t1 JR972076.1
## 4 Pocillopora_acuta_HIv2___RNAseq.g13172.t1 Gene:g13552
## 5 Pocillopora_acuta_HIv2___RNAseq.g13172.t1 aug_v2a.06327.t1
## 6 Pocillopora_acuta_HIv2___RNAseq.g13823.t1 PFX18785.1
## definition
## 1 thioredoxin reductase 1, cytoplasmic-like [Stylophora pistillata]
## 2 Flagellar associated protein
## 3 Acidic skeletal organic matrix protein (Acidic SOMP)
## 4 Acidic SOMP (Full-Length p27)
## 5 SAARP3
## 6 Mucin-4 [Stylophora pistillata]
## Ref substanceBXH
## 1 Peled et al., 2020 Pocillopora_acuta_HIv2___RNAseq.g10093.t2
## 2 Drake et al., 2013 Pocillopora_acuta_HIv2___RNAseq.g11609.t1
## 3 Ramos-Silva et al., 2013 Pocillopora_acuta_HIv2___RNAseq.g13172.t1
## 4 Mummadisetti et al., 2021 Pocillopora_acuta_HIv2___RNAseq.g13172.t1
## 5 Takeuchi et al., 2016 Pocillopora_acuta_HIv2___RNAseq.g13172.t1
## 6 Peled et al., 2020 Pocillopora_acuta_HIv2___RNAseq.g13823.t1
## geneSymbol moduleColor
## 1 Pocillopora_acuta_HIv2___Sc0000021 brown
## 2 Pocillopora_acuta_HIv2___Sc0000013 turquoise
## 3 Pocillopora_acuta_HIv2___Sc0000004 red
## 4 Pocillopora_acuta_HIv2___Sc0000004 red
## 5 Pocillopora_acuta_HIv2___Sc0000004 red
## 6 Pocillopora_acuta_HIv2___Sc0000005 pink
## GO.terms
## 1 GO:0000003,GO:0000302,GO:0000305,GO:0001650,GO:0001704,GO:0001707,GO:0001887,GO:0001890,GO:0003006,GO:0003674,GO:0003824,GO:0004791,GO:0005488,GO:0005515,GO:0005575,GO:0005622,GO:0005623,GO:0005634,GO:0005654,GO:0005730,GO:0005737,GO:0005739,GO:0005783,GO:0005829,GO:0006082,GO:0006139,GO:0006518,GO:0006520,GO:0006575,GO:0006725,GO:0006732,GO:0006733,GO:0006739,GO:0006749,GO:0006753,GO:0006790,GO:0006793,GO:0006796,GO:0006807,GO:0006950,GO:0006979,GO:0007154,GO:0007165,GO:0007275,GO:0007369,GO:0007498,GO:0008150,GO:0008152,GO:0008283,GO:0009056,GO:0009069,GO:0009117,GO:0009611,GO:0009628,GO:0009636,GO:0009653,GO:0009790,GO:0009888,GO:0009987,GO:0010035,GO:0010038,GO:0010269,GO:0010941,GO:0010942,GO:0012505,GO:0015036,GO:0015949,GO:0016043,GO:0016174,GO:0016209,GO:0016259,GO:0016491,GO:0016651,GO:0016667,GO:0016668,GO:0016999,GO:0017001,GO:0017144,GO:0018996,GO:0019216,GO:0019222,GO:0019362,GO:0019637,GO:0019725,GO:0019752,GO:0022404,GO:0022414,GO:0022607,GO:0023052,GO:0031974,GO:0031981,GO:0032501,GO:0032502,GO:0033554,GO:0033797,GO:0034599,GO:0034641,GO:0036295,GO:0036296,GO:0036477,GO:0042221,GO:0042303,GO:0042395,GO:0042493,GO:0042537,GO:0042592,GO:0042737,GO:0042743,GO:0042744,GO:0042802,GO:0042803,GO:0043025,GO:0043167,GO:0043169,GO:0043226,GO:0043227,GO:0043228,GO:0043229,GO:0043231,GO:0043232,GO:0043233,GO:0043436,GO:0043603,GO:0043933,GO:0044085,GO:0044237,GO:0044238,GO:0044248,GO:0044281,GO:0044297,GO:0044422,GO:0044424,GO:0044428,GO:0044444,GO:0044446,GO:0044452,GO:0044464,GO:0045340,GO:0045454,GO:0046483,GO:0046496,GO:0046688,GO:0046872,GO:0046914,GO:0046983,GO:0048332,GO:0048513,GO:0048518,GO:0048522,GO:0048598,GO:0048608,GO:0048646,GO:0048678,GO:0048729,GO:0048731,GO:0048856,GO:0050664,GO:0050789,GO:0050794,GO:0050896,GO:0051186,GO:0051187,GO:0051259,GO:0051262,GO:0051716,GO:0055086,GO:0055093,GO:0055114,GO:0061458,GO:0065003,GO:0065007,GO:0065008,GO:0070013,GO:0070276,GO:0070482,GO:0070887,GO:0070995,GO:0071241,GO:0071248,GO:0071280,GO:0071453,GO:0071455,GO:0071704,GO:0071840,GO:0072524,GO:0072593,GO:0080090,GO:0097237,GO:0097458,GO:0098623,GO:0098625,GO:0098626,GO:0098754,GO:0098869,GO:1901360,GO:1901564,GO:1901605,GO:1901700,GO:1990748
## 2 -
## 3 <NA>
## 4 <NA>
## 5 <NA>
## 6 <NA>
## GO.description GS.Flat GS.Slope p.GS.Flat
## 1 thioredoxin-disulfide reductase activity 0.57178848 -0.57178848 3.311055e-05
## 2 - -0.29586493 0.29586493 4.589336e-02
## 3 <NA> 0.35628512 -0.35628512 1.508700e-02
## 4 <NA> 0.35628512 -0.35628512 1.508700e-02
## 5 <NA> 0.35628512 -0.35628512 1.508700e-02
## 6 <NA> -0.05455251 0.05455251 7.187880e-01
## p.GS.Slope A.brown p.A.brown A.magenta p.A.magenta A.red
## 1 3.311055e-05 0.7005073 5.973619e-08 -0.3738439 1.048844e-02 0.2901298
## 2 4.589336e-02 -0.4291375 2.921081e-03 0.3115539 3.505853e-02 -0.3452015
## 3 1.508700e-02 0.4914202 5.241968e-04 -0.6288605 2.864308e-06 0.6673892
## 4 1.508700e-02 0.4914202 5.241968e-04 -0.6288605 2.864308e-06 0.6673892
## 5 1.508700e-02 0.4914202 5.241968e-04 -0.6288605 2.864308e-06 0.6673892
## 6 7.187880e-01 0.0972208 5.203783e-01 -0.3252127 2.743096e-02 0.3709019
## p.A.red A.turquoise p.A.turquoise A.purple p.A.purple A.green
## 1 5.047695e-02 -0.43323233 2.634293e-03 0.6984202 6.792759e-08 0.4574538
## 2 1.879503e-02 0.58815287 1.720729e-05 -0.1784560 2.353887e-01 -0.1306835
## 3 4.071016e-07 -0.13892006 3.571825e-01 0.1198762 4.274677e-01 0.2378899
## 4 4.071016e-07 -0.13892006 3.571825e-01 0.1198762 4.274677e-01 0.2378899
## 5 4.071016e-07 -0.13892006 3.571825e-01 0.1198762 4.274677e-01 0.2378899
## 6 1.116224e-02 0.08164806 5.895928e-01 -0.1391597 3.563456e-01 0.1614616
## p.A.green A.lightcyan p.A.lightcyan A.pink p.A.pink A.blue
## 1 0.001391986 -0.3508191 1.682948e-02 0.1707384 2.565893e-01 0.12358439
## 2 0.386672688 0.1196505 4.283449e-01 -0.1522331 3.125037e-01 -0.58598406
## 3 0.111386103 -0.6473842 1.159989e-06 0.7188738 1.835918e-08 0.07448551
## 4 0.111386103 -0.6473842 1.159989e-06 0.7188738 1.835918e-08 0.07448551
## 5 0.111386103 -0.6473842 1.159989e-06 0.7188738 1.835918e-08 0.07448551
## 6 0.283719167 -0.5276145 1.646022e-04 0.6417477 1.537006e-06 -0.02286640
## p.A.blue A.salmon p.A.salmon A.midnightblue p.A.midnightblue
## 1 4.132051e-01 0.1178467 0.435389343 0.2439890 0.10224333
## 2 1.880492e-05 0.1907995 0.204028320 0.2258109 0.13131383
## 3 6.227429e-01 0.4254256 0.003204458 0.2691022 0.07053914
## 4 6.227429e-01 0.4254256 0.003204458 0.2691022 0.07053914
## 5 6.227429e-01 0.4254256 0.003204458 0.2691022 0.07053914
## 6 8.801027e-01 0.2940377 0.047315397 0.2906592 0.05003895
## A.black p.A.black A.cyan p.A.cyan A.yellow p.A.yellow
## 1 -0.28430645 0.05550307 0.04904562 0.7461773 0.05522073 0.7154873547
## 2 -0.18825739 0.21023361 0.07386502 0.6256510 -0.14392338 0.3399558326
## 3 0.09618758 0.52484209 0.16699226 0.2673276 -0.38010677 0.0091694889
## 4 0.09618758 0.52484209 0.16699226 0.2673276 -0.38010677 0.0091694889
## 5 0.09618758 0.52484209 0.16699226 0.2673276 -0.38010677 0.0091694889
## 6 0.02103556 0.88964060 -0.13338389 0.3768501 -0.46998526 0.0009821983
## A.tan p.A.tan
## 1 0.2648346 0.075293267
## 2 0.2613466 0.079363532
## 3 -0.2055446 0.170565420
## 4 -0.2055446 0.170565420
## 5 -0.2055446 0.170565420
## 6 -0.3805358 0.009084622glmmseq_exp <-read.csv("../../output/glmmseq/model_expression_prediction_allgenes.csv")
glmmseq_exp <- plyr::rename(glmmseq_exp, c("X"="Pocillopora_acuta_best_hit"))
head(glmmseq_exp)## Pocillopora_acuta_best_hit y_Stable_Slope y_Variable_Slope
## 1 Pocillopora_acuta_HIv2___RNAseq.g27841.t1 47.310129 52.236750
## 2 Pocillopora_acuta_HIv2___RNAseq.g14011.t1 12.583333 10.093750
## 3 Pocillopora_acuta_HIv2___RNAseq.g7479.t3a 29.250000 27.416667
## 4 Pocillopora_acuta_HIv2___RNAseq.g7479.t3b 25.697917 21.916667
## 5 Pocillopora_acuta_HIv2___RNAseq.g3340.t1 1.770834 5.083334
## 6 Pocillopora_acuta_HIv2___RNAseq.g682.t1 39.883569 33.703575
## y_Stable_Flat y_Variable_Flat LCI_Stable_Slope LCI_Variable_Slope
## 1 55.181979 56.54215 40.0331592 44.306916
## 2 7.375000 8.37500 9.0538343 7.194628
## 3 33.454545 33.90909 25.1660941 23.531859
## 4 61.636364 53.45455 15.8852724 13.517735
## 5 8.875001 11.22727 0.4082257 1.220526
## 6 32.709968 34.68832 30.4669900 25.830257
## LCI_Stable_Flat LCI_Variable_Flat UCI_Stable_Slope UCI_Variable_Slope
## 1 46.526240 47.595803 55.909859 61.58583
## 2 5.087888 5.822403 17.488754 14.16109
## 3 28.730827 29.134602 33.996634 31.94281
## 4 37.592184 32.573563 41.572024 35.53408
## 5 2.019801 2.562261 7.681661 21.17143
## 6 24.814117 26.322279 52.210576 43.97676
## UCI_Stable_Flat UCI_Variable_Flat
## 1 65.44803 67.17011
## 2 10.69022 12.04668
## 3 38.95490 39.46601
## 4 101.05934 87.72109
## 5 38.99674 49.19549
## 6 43.11828 45.71334biomin_exp <- merge(biomin_mod, glmmseq_exp, by=c("Pocillopora_acuta_best_hit"), all=T)
head(biomin_exp)## Pocillopora_acuta_best_hit accessionnumber.geneID definition
## 1 Pocillopora_acuta_HIv2___RNAseq.10002_t <NA> <NA>
## 2 Pocillopora_acuta_HIv2___RNAseq.10171_t <NA> <NA>
## 3 Pocillopora_acuta_HIv2___RNAseq.10263_t <NA> <NA>
## 4 Pocillopora_acuta_HIv2___RNAseq.10431_t <NA> <NA>
## 5 Pocillopora_acuta_HIv2___RNAseq.10474_t <NA> <NA>
## 6 Pocillopora_acuta_HIv2___RNAseq.10518_t <NA> <NA>
## Ref substanceBXH geneSymbol moduleColor GO.terms GO.description GS.Flat
## 1 <NA> <NA> <NA> <NA> <NA> <NA> NA
## 2 <NA> <NA> <NA> <NA> <NA> <NA> NA
## 3 <NA> <NA> <NA> <NA> <NA> <NA> NA
## 4 <NA> <NA> <NA> <NA> <NA> <NA> NA
## 5 <NA> <NA> <NA> <NA> <NA> <NA> NA
## 6 <NA> <NA> <NA> <NA> <NA> <NA> NA
## GS.Slope p.GS.Flat p.GS.Slope A.brown p.A.brown A.magenta p.A.magenta A.red
## 1 NA NA NA NA NA NA NA NA
## 2 NA NA NA NA NA NA NA NA
## 3 NA NA NA NA NA NA NA NA
## 4 NA NA NA NA NA NA NA NA
## 5 NA NA NA NA NA NA NA NA
## 6 NA NA NA NA NA NA NA NA
## p.A.red A.turquoise p.A.turquoise A.purple p.A.purple A.green p.A.green
## 1 NA NA NA NA NA NA NA
## 2 NA NA NA NA NA NA NA
## 3 NA NA NA NA NA NA NA
## 4 NA NA NA NA NA NA NA
## 5 NA NA NA NA NA NA NA
## 6 NA NA NA NA NA NA NA
## A.lightcyan p.A.lightcyan A.pink p.A.pink A.blue p.A.blue A.salmon p.A.salmon
## 1 NA NA NA NA NA NA NA NA
## 2 NA NA NA NA NA NA NA NA
## 3 NA NA NA NA NA NA NA NA
## 4 NA NA NA NA NA NA NA NA
## 5 NA NA NA NA NA NA NA NA
## 6 NA NA NA NA NA NA NA NA
## A.midnightblue p.A.midnightblue A.black p.A.black A.cyan p.A.cyan A.yellow
## 1 NA NA NA NA NA NA NA
## 2 NA NA NA NA NA NA NA
## 3 NA NA NA NA NA NA NA
## 4 NA NA NA NA NA NA NA
## 5 NA NA NA NA NA NA NA
## 6 NA NA NA NA NA NA NA
## p.A.yellow A.tan p.A.tan y_Stable_Slope y_Variable_Slope y_Stable_Flat
## 1 NA NA NA 9.291667 7.302083 9.295455
## 2 NA NA NA 11.718750 18.281250 116.204506
## 3 NA NA NA 342.041619 416.236831 460.406057
## 4 NA NA NA 13.468750 12.614583 10.318182
## 5 NA NA NA 45.079954 39.602336 23.349020
## 6 NA NA NA 24.020833 16.291667 18.863636
## y_Variable_Flat LCI_Stable_Slope LCI_Variable_Slope LCI_Stable_Flat
## 1 9.556818 4.926111 3.842939 4.791061
## 2 31.727281 5.165294 8.106281 50.173956
## 3 502.653747 288.140041 351.522240 386.585819
## 4 14.295455 6.276418 5.872199 4.624811
## 5 24.533092 25.255593 22.561811 9.268724
## 6 22.204545 11.549652 7.799681 8.756677
## LCI_Variable_Flat UCI_Stable_Slope UCI_Variable_Slope UCI_Stable_Flat
## 1 4.929596 17.52601 13.87491 18.03473
## 2 13.633889 26.58689 41.22780 269.13340
## 3 421.462603 406.02642 492.86526 548.32259
## 4 6.445353 28.90298 27.09849 23.02038
## 5 9.999326 80.46543 69.51326 58.81896
## 6 10.326296 49.95825 34.02939 40.63605
## UCI_Variable_Flat
## 1 18.52744
## 2 73.83222
## 3 599.48567
## 4 31.70657
## 5 60.19132
## 6 47.74624biomin_all <- read.csv("../../output/Biomin_blast_Pocillopora_acuta_best_hit_glmmSeq.csv")
head(biomin_all)## Gene Dispersion AIC logLik
## 1 Pocillopora_acuta_HIv2___RNAseq.g13823.t1 0.29270610 366.4332 -177.2166
## 2 Pocillopora_acuta_HIv2___RNAseq.g13823.t1 0.29270610 366.4332 -177.2166
## 3 Pocillopora_acuta_HIv2___RNAseq.g13823.t1 0.29270610 366.4332 -177.2166
## 4 Pocillopora_acuta_HIv2___RNAseq.g25351.t1 0.13435721 815.9720 -401.9860
## 5 Pocillopora_acuta_HIv2___RNAseq.g7085.t1 0.25162583 491.9378 -239.9689
## 6 Pocillopora_acuta_HIv2___RNAseq.g22851.t1 0.02235176 759.7088 -373.8544
## meanExp X.Intercept. TreatmentVariable OriginFlat
## 1 4.316132 3.079614 0.16922532 0.07545055
## 2 4.316132 3.079614 0.16922532 0.07545055
## 3 4.316132 3.079614 0.16922532 0.07545055
## 4 12.237463 8.288576 -0.01709198 0.40440178
## 5 6.477549 4.543751 -0.19119978 0.09230384
## 6 11.722317 8.185091 -0.11783143 -0.07240262
## TreatmentVariable.OriginFlat se_.Intercept. se_TreatmentVariable
## 1 -0.3792527 0.1679989 0.2363277
## 2 -0.3792527 0.1679989 0.2363277
## 3 -0.3792527 0.1679989 0.2363277
## 4 0.2217089 0.1059042 0.1497752
## 5 0.3564220 0.1641296 0.2123867
## 6 0.1972383 0.0571238 0.0618745
## se_OriginFlat se_TreatmentVariable.OriginFlat Chisq_Treatment Chisq_Origin
## 1 0.24229985 0.34311487 0.003896091 0.4390691
## 2 0.24229985 0.34311487 0.003896091 0.4390691
## 3 0.24229985 0.34311487 0.003896091 0.4390691
## 4 0.15311285 0.21653465 0.676680451 22.6483669
## 5 0.22652466 0.30543197 0.013908932 2.6586173
## 6 0.08256803 0.08945199 0.275523829 0.1398413
## Chisq_Treatment.Origin Df_Treatment Df_Origin Df_Treatment.Origin P_Treatment
## 1 1.221739 1 1 1 0.9502294
## 2 1.221739 1 1 1 0.9502294
## 3 1.221739 1 1 1 0.9502294
## 4 1.048362 1 1 1 0.4107321
## 5 1.361758 1 1 1 0.9061183
## 6 4.861862 1 1 1 0.5996502
## P_Origin P_Treatment.Origin Treatment Origin Treatment.Origin
## 1 5.075721e-01 0.26901967 1 0.5306740470 0.9994279
## 2 5.075721e-01 0.26901967 1 0.5306740470 0.9994279
## 3 5.075721e-01 0.26901967 1 0.5306740470 0.9994279
## 4 1.945255e-06 0.30588455 1 0.0001683701 0.9994279
## 5 1.029902e-01 0.24323297 1 0.2466098058 0.9994279
## 6 7.084388e-01 0.02745669 1 0.6038776397 0.9994279
## Stable_OriginFC Variable_OriginFC maxGroupFC col RF13B
## 1 -0.1042361 0.4192811 Variable Not Significant 21.57376
## 2 -0.1042361 0.4192811 Variable Not Significant 21.57376
## 3 -0.1042361 0.4192811 Variable Not Significant 21.57376
## 4 -0.5833078 -0.9031151 Variable q_Origin < 0.05 5832.10758
## 5 -0.1318273 -0.6407294 Variable Not Significant 67.80326
## 6 0.1044247 -0.1800468 Variable Not Significant 4739.03686
## RF13D RF14B RF14C RF15B RF15D RF17B RF17D
## 1 25.50848 16.94826 19.20733 28.23482 21.29134 25.87585 14.01704
## 2 25.50848 16.94826 19.20733 28.23482 21.29134 25.87585 14.01704
## 3 25.50848 16.94826 19.20733 28.23482 21.29134 25.87585 14.01704
## 4 4927.06132 13021.91653 6113.26692 5054.03248 7484.79414 8321.89749 7240.96832
## 5 70.63887 54.61107 99.23788 164.36698 112.66669 121.50398 91.11075
## 6 4597.41325 2757.85926 3654.72843 3468.84911 3098.77752 2764.21551 2954.09080
## RF18B RF18D RF19B RF19C RF20B RF20C RF22B
## 1 22.30905 43.47845 0.0000 24.20314 34.69447 16.42775 13.92415
## 2 22.30905 43.47845 0.0000 24.20314 34.69447 16.42775 13.92415
## 3 22.30905 43.47845 0.0000 24.20314 34.69447 16.42775 13.92415
## 4 5891.93841 5089.20850 8016.9416 8866.41520 4982.12595 8081.53904 7060.31604
## 5 123.28687 101.44972 118.3772 140.60869 115.64823 135.98524 42.54600
## 6 2417.59689 3150.51549 3524.0588 2869.80032 3728.49908 2716.05423 4221.33720
## RF22C RF23A RF23C RF24B RF24D RF25A RF25C
## 1 18.19144 18.60133 36.2419 18.57332 23.35614 10.69681 16.36654
## 2 18.19144 18.60133 36.2419 18.57332 23.35614 10.69681 16.36654
## 3 18.19144 18.60133 36.2419 18.57332 23.35614 10.69681 16.36654
## 4 3512.96929 6530.04515 4346.0891 4820.26566 5124.51037 6260.69212 4758.11817
## 5 46.48924 249.64940 184.1480 184.75567 70.06842 92.45103 98.19927
## 6 3153.18302 3553.83267 3132.4753 3154.53200 3760.33872 3892.11199 3634.28212
## RS11B RS11D RS12A RS12C RS13A RS13C RS14B
## 1 37.03275 14.09808 0.00000 8.755265 20.55198 15.25056 29.16535
## 2 37.03275 14.09808 0.00000 8.755265 20.55198 15.25056 29.16535
## 3 37.03275 14.09808 0.00000 8.755265 20.55198 15.25056 29.16535
## 4 4408.84668 3804.31239 2168.53312 3681.588844 2701.11675 2417.21364 3177.56482
## 5 75.04005 300.39753 20.40972 12.257371 104.71721 80.06544 85.30865
## 6 3502.51878 3288.10581 3286.98596 2969.785817 3332.35599 3746.87177 3115.58845
## RS14C RS15B RS15D RS1B RS1C RS2B RS2C
## 1 14.95038 17.20235 15.77586 6.823476 33.59846 34.30448 35.76535
## 2 14.95038 17.20235 15.77586 6.823476 33.59846 34.30448 35.76535
## 3 14.95038 17.20235 15.77586 6.823476 33.59846 34.30448 35.76535
## 4 3123.75088 2933.85998 3892.69273 3421.973109 3717.90918 6782.61209 5170.97795
## 5 73.87249 85.15161 46.01292 55.725052 95.99559 59.81295 85.37536
## 6 3313.70871 3593.56992 3662.62815 3472.011931 3480.80008 2515.66212 3603.07098
## RS3B RS3D RS6A RS6D RS7B RS7C RS8B
## 1 44.72649 29.74233 29.32862 38.32429 30.78679 17.29086 9.894721
## 2 44.72649 29.74233 29.32862 38.32429 30.78679 17.29086 9.894721
## 3 44.72649 29.74233 29.32862 38.32429 30.78679 17.29086 9.894721
## 4 5164.93773 4277.13888 5058.41471 4024.05055 2932.44205 5029.82092 3791.477100
## 5 118.62244 77.71383 61.74446 120.29569 71.46934 82.81413 52.172164
## 6 2804.15671 2625.00044 3068.69956 3619.51637 3201.82649 3680.22360 3224.779455
## RS8C RS9A RS9C accessionnumber.geneID
## 1 14.95208 47.19477 19.44197 aug_v2a.09809.t1
## 2 14.95208 47.19477 19.44197 P13_g6918
## 3 14.95208 47.19477 19.44197 PFX18785.1
## 4 4200.60065 4097.88734 3391.45658 XP_022794351.1
## 5 75.69492 145.03759 102.65358 XP_022799541.1
## 6 3453.93103 2906.50718 4526.08973 P4_g9861
## definition
## 1 Mucin4-like protein
## 2 Sushi domain-containing
## 3 Mucin-4 [Stylophora pistillata]
## 4 mammalian ependymin-related protein 1-like [Stylophora pistillata]
## 5 uncharacterized protein LOC111337489 [Stylophora pistillata]
## 6 Viral inclusion protein
## Ref
## 1 Takeuchi et al., 2016
## 2 Drake et al., 2013
## 3 Peled et al., 2020
## 4 Peled et al., 2020
## 5 Peled et al., 2020
## 6 Drake et al., 2013colnames(biomin_all)## [1] "Gene" "Dispersion"
## [3] "AIC" "logLik"
## [5] "meanExp" "X.Intercept."
## [7] "TreatmentVariable" "OriginFlat"
## [9] "TreatmentVariable.OriginFlat" "se_.Intercept."
## [11] "se_TreatmentVariable" "se_OriginFlat"
## [13] "se_TreatmentVariable.OriginFlat" "Chisq_Treatment"
## [15] "Chisq_Origin" "Chisq_Treatment.Origin"
## [17] "Df_Treatment" "Df_Origin"
## [19] "Df_Treatment.Origin" "P_Treatment"
## [21] "P_Origin" "P_Treatment.Origin"
## [23] "Treatment" "Origin"
## [25] "Treatment.Origin" "Stable_OriginFC"
## [27] "Variable_OriginFC" "maxGroupFC"
## [29] "col" "RF13B"
## [31] "RF13D" "RF14B"
## [33] "RF14C" "RF15B"
## [35] "RF15D" "RF17B"
## [37] "RF17D" "RF18B"
## [39] "RF18D" "RF19B"
## [41] "RF19C" "RF20B"
## [43] "RF20C" "RF22B"
## [45] "RF22C" "RF23A"
## [47] "RF23C" "RF24B"
## [49] "RF24D" "RF25A"
## [51] "RF25C" "RS11B"
## [53] "RS11D" "RS12A"
## [55] "RS12C" "RS13A"
## [57] "RS13C" "RS14B"
## [59] "RS14C" "RS15B"
## [61] "RS15D" "RS1B"
## [63] "RS1C" "RS2B"
## [65] "RS2C" "RS3B"
## [67] "RS3D" "RS6A"
## [69] "RS6D" "RS7B"
## [71] "RS7C" "RS8B"
## [73] "RS8C" "RS9A"
## [75] "RS9C" "accessionnumber.geneID"
## [77] "definition" "Ref"library(tidyr)
biomin_all_filtered_long <- pivot_longer(biomin_all, cols=30:75, names_to = "Colony", values_to = "Expression")
biomin_all_filtered_long $Colony <- as.factor(biomin_all_filtered_long $Colony)
head(biomin_all_filtered_long)## # A tibble: 6 × 34
## Gene Dispersion AIC logLik meanExp X.Intercept. TreatmentVariable
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Pocillopora_ac… 0.293 366. -177. 4.32 3.08 0.169
## 2 Pocillopora_ac… 0.293 366. -177. 4.32 3.08 0.169
## 3 Pocillopora_ac… 0.293 366. -177. 4.32 3.08 0.169
## 4 Pocillopora_ac… 0.293 366. -177. 4.32 3.08 0.169
## 5 Pocillopora_ac… 0.293 366. -177. 4.32 3.08 0.169
## 6 Pocillopora_ac… 0.293 366. -177. 4.32 3.08 0.169
## # ℹ 27 more variables: OriginFlat <dbl>, TreatmentVariable.OriginFlat <dbl>,
## # se_.Intercept. <dbl>, se_TreatmentVariable <dbl>, se_OriginFlat <dbl>,
## # se_TreatmentVariable.OriginFlat <dbl>, Chisq_Treatment <dbl>,
## # Chisq_Origin <dbl>, Chisq_Treatment.Origin <dbl>, Df_Treatment <int>,
## # Df_Origin <int>, Df_Treatment.Origin <int>, P_Treatment <dbl>,
## # P_Origin <dbl>, P_Treatment.Origin <dbl>, Treatment <dbl>, Origin <dbl>,
## # Treatment.Origin <dbl>, Stable_OriginFC <dbl>, Variable_OriginFC <dbl>, …biomin_all_filtered_long <- biomin_all_filtered_long %>%
separate(Colony, into = c('Origin', 'Colony.number'), sep = 2)
head(biomin_all_filtered_long)## # A tibble: 6 × 34
## Gene Dispersion AIC logLik meanExp X.Intercept. TreatmentVariable
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Pocillopora_ac… 0.293 366. -177. 4.32 3.08 0.169
## 2 Pocillopora_ac… 0.293 366. -177. 4.32 3.08 0.169
## 3 Pocillopora_ac… 0.293 366. -177. 4.32 3.08 0.169
## 4 Pocillopora_ac… 0.293 366. -177. 4.32 3.08 0.169
## 5 Pocillopora_ac… 0.293 366. -177. 4.32 3.08 0.169
## 6 Pocillopora_ac… 0.293 366. -177. 4.32 3.08 0.169
## # ℹ 27 more variables: OriginFlat <dbl>, TreatmentVariable.OriginFlat <dbl>,
## # se_.Intercept. <dbl>, se_TreatmentVariable <dbl>, se_OriginFlat <dbl>,
## # se_TreatmentVariable.OriginFlat <dbl>, Chisq_Treatment <dbl>,
## # Chisq_Origin <dbl>, Chisq_Treatment.Origin <dbl>, Df_Treatment <int>,
## # Df_Origin <int>, Df_Treatment.Origin <int>, P_Treatment <dbl>,
## # P_Origin <dbl>, P_Treatment.Origin <dbl>, Treatment <dbl>,
## # Treatment.Origin <dbl>, Stable_OriginFC <dbl>, Variable_OriginFC <dbl>, …library(stringr)
biomin_all_filtered_long $Colony <- as.numeric(str_extract(biomin_all_filtered_long $Colony.number, "[0-9]+"))
biomin_all_filtered_long <-biomin_all_filtered_long %>%
mutate(Treatment = trimws(str_remove(biomin_all_filtered_long $Colony.number, "(\\s+[A-Za-z]+)?[0-9-]+")))
head(biomin_all_filtered_long)## # A tibble: 6 × 35
## Gene Dispersion AIC logLik meanExp X.Intercept. TreatmentVariable
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Pocillopora_ac… 0.293 366. -177. 4.32 3.08 0.169
## 2 Pocillopora_ac… 0.293 366. -177. 4.32 3.08 0.169
## 3 Pocillopora_ac… 0.293 366. -177. 4.32 3.08 0.169
## 4 Pocillopora_ac… 0.293 366. -177. 4.32 3.08 0.169
## 5 Pocillopora_ac… 0.293 366. -177. 4.32 3.08 0.169
## 6 Pocillopora_ac… 0.293 366. -177. 4.32 3.08 0.169
## # ℹ 28 more variables: OriginFlat <dbl>, TreatmentVariable.OriginFlat <dbl>,
## # se_.Intercept. <dbl>, se_TreatmentVariable <dbl>, se_OriginFlat <dbl>,
## # se_TreatmentVariable.OriginFlat <dbl>, Chisq_Treatment <dbl>,
## # Chisq_Origin <dbl>, Chisq_Treatment.Origin <dbl>, Df_Treatment <int>,
## # Df_Origin <int>, Df_Treatment.Origin <int>, P_Treatment <dbl>,
## # P_Origin <dbl>, P_Treatment.Origin <dbl>, Treatment <chr>,
## # Treatment.Origin <dbl>, Stable_OriginFC <dbl>, Variable_OriginFC <dbl>, …biomin_all_filtered_long $Origin <- as.factor(biomin_all_filtered_long $Origin)
biomin_all_filtered_long $Treatment <- as.factor(biomin_all_filtered_long $Treatment)
head(biomin_all_filtered_long)## # A tibble: 6 × 35
## Gene Dispersion AIC logLik meanExp X.Intercept. TreatmentVariable
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Pocillopora_ac… 0.293 366. -177. 4.32 3.08 0.169
## 2 Pocillopora_ac… 0.293 366. -177. 4.32 3.08 0.169
## 3 Pocillopora_ac… 0.293 366. -177. 4.32 3.08 0.169
## 4 Pocillopora_ac… 0.293 366. -177. 4.32 3.08 0.169
## 5 Pocillopora_ac… 0.293 366. -177. 4.32 3.08 0.169
## 6 Pocillopora_ac… 0.293 366. -177. 4.32 3.08 0.169
## # ℹ 28 more variables: OriginFlat <dbl>, TreatmentVariable.OriginFlat <dbl>,
## # se_.Intercept. <dbl>, se_TreatmentVariable <dbl>, se_OriginFlat <dbl>,
## # se_TreatmentVariable.OriginFlat <dbl>, Chisq_Treatment <dbl>,
## # Chisq_Origin <dbl>, Chisq_Treatment.Origin <dbl>, Df_Treatment <int>,
## # Df_Origin <int>, Df_Treatment.Origin <int>, P_Treatment <dbl>,
## # P_Origin <dbl>, P_Treatment.Origin <dbl>, Treatment <fct>,
## # Treatment.Origin <dbl>, Stable_OriginFC <dbl>, Variable_OriginFC <dbl>, …biomin_all_filtered_long <- biomin_all_filtered_long %>%
mutate(Treatment2 = ifelse(Treatment == "A" | Treatment == "B", "Variable",
ifelse(Treatment == "C" | Treatment == "D", "Stable", NA)))
biomin_all_filtered_long $Treatment2 <- as.factor(biomin_all_filtered_long $Treatment2)
head(biomin_all_filtered_long)## # A tibble: 6 × 36
## Gene Dispersion AIC logLik meanExp X.Intercept. TreatmentVariable
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Pocillopora_ac… 0.293 366. -177. 4.32 3.08 0.169
## 2 Pocillopora_ac… 0.293 366. -177. 4.32 3.08 0.169
## 3 Pocillopora_ac… 0.293 366. -177. 4.32 3.08 0.169
## 4 Pocillopora_ac… 0.293 366. -177. 4.32 3.08 0.169
## 5 Pocillopora_ac… 0.293 366. -177. 4.32 3.08 0.169
## 6 Pocillopora_ac… 0.293 366. -177. 4.32 3.08 0.169
## # ℹ 29 more variables: OriginFlat <dbl>, TreatmentVariable.OriginFlat <dbl>,
## # se_.Intercept. <dbl>, se_TreatmentVariable <dbl>, se_OriginFlat <dbl>,
## # se_TreatmentVariable.OriginFlat <dbl>, Chisq_Treatment <dbl>,
## # Chisq_Origin <dbl>, Chisq_Treatment.Origin <dbl>, Df_Treatment <int>,
## # Df_Origin <int>, Df_Treatment.Origin <int>, P_Treatment <dbl>,
## # P_Origin <dbl>, P_Treatment.Origin <dbl>, Treatment <fct>,
## # Treatment.Origin <dbl>, Stable_OriginFC <dbl>, Variable_OriginFC <dbl>, …MERed signficant genes
red_genes <- biomin_mod %>% filter(moduleColor == "red") %>% pull(Pocillopora_acuta_best_hit) %>% unique()
length(red_genes)## [1] 7red_genes_data <- biomin_all %>% filter(Gene %in% red_genes)
red_genes_data## Gene Dispersion AIC logLik
## 1 Pocillopora_acuta_HIv2___RNAseq.g15280.t1 0.3761588 424.9745 -206.4872
## 2 Pocillopora_acuta_HIv2___RNAseq.g13172.t1 2.2644800 299.6692 -143.8346
## 3 Pocillopora_acuta_HIv2___RNAseq.g13172.t1 2.2644800 299.6692 -143.8346
## 4 Pocillopora_acuta_HIv2___RNAseq.g13172.t1 2.2644800 299.6692 -143.8346
## 5 Pocillopora_acuta_HIv2___RNAseq.g26037.t1 1.0123067 357.7130 -172.8565
## 6 Pocillopora_acuta_HIv2___RNAseq.g7402.t1 0.2289670 640.4572 -314.2286
## 7 Pocillopora_acuta_HIv2___TS.g12304.t1 0.4572350 841.4733 -414.7366
## 8 Pocillopora_acuta_HIv2___TS.g12304.t1 0.4572350 841.4733 -414.7366
## 9 Pocillopora_acuta_HIv2___TS.g12304.t1 0.4572350 841.4733 -414.7366
## 10 Pocillopora_acuta_HIv2___TS.g12304.t1 0.4572350 841.4733 -414.7366
## 11 Pocillopora_acuta_HIv2___RNAseq.g7908.t1 0.1379748 426.1346 -207.0673
## 12 Pocillopora_acuta_HIv2___RNAseq.g6446.t1 0.2837901 526.0525 -257.0262
## 13 Pocillopora_acuta_HIv2___RNAseq.g6446.t1 0.2837901 526.0525 -257.0262
## 14 Pocillopora_acuta_HIv2___RNAseq.g6446.t1 0.2837901 526.0525 -257.0262
## 15 Pocillopora_acuta_HIv2___RNAseq.g6446.t1 0.2837901 526.0525 -257.0262
## 16 Pocillopora_acuta_HIv2___RNAseq.g6446.t1 0.2837901 526.0525 -257.0262
## 17 Pocillopora_acuta_HIv2___RNAseq.g6446.t1 0.2837901 526.0525 -257.0262
## 18 Pocillopora_acuta_HIv2___RNAseq.g6446.t1 0.2837901 526.0525 -257.0262
## meanExp X.Intercept. TreatmentVariable OriginFlat
## 1 5.119627 3.445015 0.08897577 0.5738221
## 2 2.418942 2.007934 -0.97016372 0.6575562
## 3 2.418942 2.007934 -0.97016372 0.6575562
## 4 2.418942 2.007934 -0.97016372 0.6575562
## 5 3.419645 2.773240 -0.22444560 0.1627910
## 6 8.970229 6.139525 0.01285273 0.4262517
## 7 11.599955 7.883666 0.04482990 0.8109426
## 8 11.599955 7.883666 0.04482990 0.8109426
## 9 11.599955 7.883666 0.04482990 0.8109426
## 10 11.599955 7.883666 0.04482990 0.8109426
## 11 5.951789 3.979682 -0.07434764 0.4391590
## 12 7.032579 4.747826 0.09109743 0.5344539
## 13 7.032579 4.747826 0.09109743 0.5344539
## 14 7.032579 4.747826 0.09109743 0.5344539
## 15 7.032579 4.747826 0.09109743 0.5344539
## 16 7.032579 4.747826 0.09109743 0.5344539
## 17 7.032579 4.747826 0.09109743 0.5344539
## 18 7.032579 4.747826 0.09109743 0.5344539
## TreatmentVariable.OriginFlat se_.Intercept. se_TreatmentVariable
## 1 -0.23661177 0.1844049 0.2603537
## 2 0.81574950 0.4470973 0.6466257
## 3 0.81574950 0.4470973 0.6466257
## 4 0.81574950 0.4470973 0.6466257
## 5 -0.05882322 0.2992723 0.4247794
## 6 -0.16546888 0.1387782 0.1962568
## 7 -0.27702814 0.1952245 0.2760692
## 8 -0.27702814 0.1952245 0.2760692
## 9 -0.27702814 0.1952245 0.2760692
## 10 -0.27702814 0.1952245 0.2760692
## 11 0.07105636 0.1142605 0.1619602
## 12 -0.28522045 0.1561142 0.2206354
## 13 -0.28522045 0.1561142 0.2206354
## 14 -0.28522045 0.1561142 0.2206354
## 15 -0.28522045 0.1561142 0.2206354
## 16 -0.28522045 0.1561142 0.2206354
## 17 -0.28522045 0.1561142 0.2206354
## 18 -0.28522045 0.1561142 0.2206354
## se_OriginFlat se_TreatmentVariable.OriginFlat Chisq_Treatment Chisq_Origin
## 1 0.2642637 0.3737703 0.01911573 5.9417187
## 2 0.6419358 0.9184510 1.51822793 5.2911915
## 3 0.6419358 0.9184510 1.51822793 5.2911915
## 4 0.6419358 0.9184510 1.51822793 5.2911915
## 5 0.4317601 0.6129627 0.68086838 0.1900494
## 6 0.2005027 0.2835874 0.21963546 5.8698042
## 7 0.2822257 0.3991048 0.19323559 11.3479034
## 8 0.2822257 0.3991048 0.19323559 11.3479034
## 9 0.2822257 0.3991048 0.19323559 11.3479034
## 10 0.2822257 0.3991048 0.19323559 11.3479034
## 11 0.1633824 0.2313253 0.11676889 16.8375086
## 12 0.2250165 0.3182766 0.08355311 6.0643147
## 13 0.2250165 0.3182766 0.08355311 6.0643147
## 14 0.2250165 0.3182766 0.08355311 6.0643147
## 15 0.2250165 0.3182766 0.08355311 6.0643147
## 16 0.2250165 0.3182766 0.08355311 6.0643147
## 17 0.2250165 0.3182766 0.08355311 6.0643147
## 18 0.2250165 0.3182766 0.08355311 6.0643147
## Chisq_Treatment.Origin Df_Treatment Df_Origin Df_Treatment.Origin
## 1 0.400740419 1 1 1
## 2 0.788863129 1 1 1
## 3 0.788863129 1 1 1
## 4 0.788863129 1 1 1
## 5 0.009209363 1 1 1
## 6 0.340454315 1 1 1
## 7 0.481807896 1 1 1
## 8 0.481807896 1 1 1
## 9 0.481807896 1 1 1
## 10 0.481807896 1 1 1
## 11 0.094353823 1 1 1
## 12 0.803067041 1 1 1
## 13 0.803067041 1 1 1
## 14 0.803067041 1 1 1
## 15 0.803067041 1 1 1
## 16 0.803067041 1 1 1
## 17 0.803067041 1 1 1
## 18 0.803067041 1 1 1
## P_Treatment P_Origin P_Treatment.Origin Treatment Origin
## 1 0.8900352 1.478659e-02 0.5267071 1 0.07884531
## 2 0.2178879 2.143355e-02 0.3744441 1 0.09820690
## 3 0.2178879 2.143355e-02 0.3744441 1 0.09820690
## 4 0.2178879 2.143355e-02 0.3744441 1 0.09820690
## 5 0.4092879 6.628755e-01 0.9235480 1 0.58754533
## 6 0.6393178 1.540277e-02 0.5595671 1 0.08113362
## 7 0.6602372 7.553319e-04 0.4876046 1 0.01164322
## 8 0.6602372 7.553319e-04 0.4876046 1 0.01164322
## 9 0.6602372 7.553319e-04 0.4876046 1 0.01164322
## 10 0.6602372 7.553319e-04 0.4876046 1 0.01164322
## 11 0.7325657 4.072046e-05 0.7587135 1 0.00158533
## 12 0.7725389 1.379403e-02 0.3701780 1 0.07587599
## 13 0.7725389 1.379403e-02 0.3701780 1 0.07587599
## 14 0.7725389 1.379403e-02 0.3701780 1 0.07587599
## 15 0.7725389 1.379403e-02 0.3701780 1 0.07587599
## 16 0.7725389 1.379403e-02 0.3701780 1 0.07587599
## 17 0.7725389 1.379403e-02 0.3701780 1 0.07587599
## 18 0.7725389 1.379403e-02 0.3701780 1 0.07587599
## Treatment.Origin Stable_OriginFC Variable_OriginFC maxGroupFC
## 1 0.9994279 -0.8082417 -0.4747322 Stable
## 2 0.9994279 -0.8639189 -1.8006515 Variable
## 3 0.9994279 -0.8639189 -1.8006515 Variable
## 4 0.9994279 -0.8639189 -1.8006515 Variable
## 5 0.9994279 -0.2220597 -0.1396272 Stable
## 6 0.9994279 -0.6138737 -0.3755265 Stable
## 7 0.9994279 -1.1696409 -0.7700607 Stable
## 8 0.9994279 -1.1696409 -0.7700607 Stable
## 9 0.9994279 -1.1696409 -0.7700607 Stable
## 10 0.9994279 -1.1696409 -0.7700607 Stable
## 11 0.9994279 -0.6241331 -0.7246607 Variable
## 12 0.9994279 -0.7659105 -0.3570662 Stable
## 13 0.9994279 -0.7659105 -0.3570662 Stable
## 14 0.9994279 -0.7659105 -0.3570662 Stable
## 15 0.9994279 -0.7659105 -0.3570662 Stable
## 16 0.9994279 -0.7659105 -0.3570662 Stable
## 17 0.9994279 -0.7659105 -0.3570662 Stable
## 18 0.9994279 -0.7659105 -0.3570662 Stable
## col RF13B RF13D RF14B RF14C RF15B
## 1 Not Significant 30.81966 33.35725 27.305537 32.012220 94.78832
## 2 Not Significant 0.00000 10.79205 2.824711 9.603666 65.54511
## 3 Not Significant 0.00000 10.79205 2.824711 9.603666 65.54511
## 4 Not Significant 0.00000 10.79205 2.824711 9.603666 65.54511
## 5 Not Significant 14.38251 44.14929 0.000000 9.603666 34.28514
## 6 Not Significant 401.68294 770.15992 308.835041 557.012625 805.70071
## 7 q_Origin < 0.05 2931.97728 4547.37738 2776.690661 4994.973370 7768.60859
## 8 q_Origin < 0.05 2931.97728 4547.37738 2776.690661 4994.973370 7768.60859
## 9 q_Origin < 0.05 2931.97728 4547.37738 2776.690661 4994.973370 7768.60859
## 10 q_Origin < 0.05 2931.97728 4547.37738 2776.690661 4994.973370 7768.60859
## 11 q_Origin < 0.05 52.39343 71.61997 95.098595 84.298846 107.89734
## 12 Not Significant 89.37702 182.48375 84.741322 153.658655 250.07982
## 13 Not Significant 89.37702 182.48375 84.741322 153.658655 250.07982
## 14 Not Significant 89.37702 182.48375 84.741322 153.658655 250.07982
## 15 Not Significant 89.37702 182.48375 84.741322 153.658655 250.07982
## 16 Not Significant 89.37702 182.48375 84.741322 153.658655 250.07982
## 17 Not Significant 89.37702 182.48375 84.741322 153.658655 250.07982
## 18 Not Significant 89.37702 182.48375 84.741322 153.658655 250.07982
## RF15D RF17B RF17D RF18B RF18D RF19B
## 1 67.42259 55.126805 45.55537 75.14629 71.34925 20.248735
## 2 17.74279 5.625184 14.01704 31.70234 30.10047 0.000000
## 3 17.74279 5.625184 14.01704 31.70234 30.10047 0.000000
## 4 17.74279 5.625184 14.01704 31.70234 30.10047 0.000000
## 5 12.41995 12.375405 18.68938 28.17986 23.41147 7.787975
## 6 765.60120 626.645519 622.59011 1063.78960 911.93265 262.454763
## 7 7687.06189 4441.645437 4540.35229 8097.01221 7354.54732 1811.483020
## 8 7687.06189 4441.645437 4540.35229 8097.01221 7354.54732 1811.483020
## 9 7687.06189 4441.645437 4540.35229 8097.01221 7354.54732 1811.483020
## 10 7687.06189 4441.645437 4540.35229 8097.01221 7354.54732 1811.483020
## 11 84.27823 68.627247 87.60649 108.02278 75.80858 76.322156
## 12 224.44624 146.254789 165.86829 290.01769 219.62192 93.455702
## 13 224.44624 146.254789 165.86829 290.01769 219.62192 93.455702
## 14 224.44624 146.254789 165.86829 290.01769 219.62192 93.455702
## 15 224.44624 146.254789 165.86829 290.01769 219.62192 93.455702
## 16 224.44624 146.254789 165.86829 290.01769 219.62192 93.455702
## 17 224.44624 146.254789 165.86829 290.01769 219.62192 93.455702
## 18 224.44624 146.254789 165.86829 290.01769 219.62192 93.455702
## RF19C RF20B RF20C RF22B RF22C RF23A
## 1 72.60941 72.85839 31.03019 26.301166 100.05292 40.139708
## 2 10.37277 18.50372 10.95183 0.000000 18.19144 7.832138
## 3 10.37277 18.50372 10.95183 0.000000 18.19144 7.832138
## 4 10.37277 18.50372 10.95183 0.000000 18.19144 7.832138
## 5 0.00000 10.40834 15.51509 9.282765 18.19144 11.748207
## 6 524.40126 1237.43611 513.82343 361.254256 1105.63533 487.550597
## 7 5988.54717 9398.73203 4582.42885 3287.645797 7853.64912 4238.165735
## 8 5988.54717 9398.73203 4582.42885 3287.645797 7853.64912 4238.165735
## 9 5988.54717 9398.73203 4582.42885 3287.645797 7853.64912 4238.165735
## 10 5988.54717 9398.73203 4582.42885 3287.645797 7853.64912 4238.165735
## 11 88.74483 83.26673 71.18691 80.450627 117.23373 97.901726
## 12 191.32002 326.12802 134.15994 99.016156 259.73334 176.223108
## 13 191.32002 326.12802 134.15994 99.016156 259.73334 176.223108
## 14 191.32002 326.12802 134.15994 99.016156 259.73334 176.223108
## 15 191.32002 326.12802 134.15994 99.016156 259.73334 176.223108
## 16 191.32002 326.12802 134.15994 99.016156 259.73334 176.223108
## 17 191.32002 326.12802 134.15994 99.016156 259.73334 176.223108
## 18 191.32002 326.12802 134.15994 99.016156 259.73334 176.223108
## RF23C RF24B RF24D RF25A RF25C RS11B
## 1 85.21743 58.65259 49.30741 29.034208 26.36832 46.778214
## 2 35.26239 10.75298 0.00000 0.000000 3.63701 7.796369
## 3 35.26239 10.75298 0.00000 0.000000 3.63701 7.796369
## 4 35.26239 10.75298 0.00000 0.000000 3.63701 7.796369
## 5 44.07798 23.46104 0.00000 4.584349 20.00355 24.363653
## 6 1178.35140 642.24590 519.02536 533.312558 343.69743 547.694924
## 7 11209.52075 4778.23130 5289.73344 2901.128629 1513.90536 5039.378029
## 8 11209.52075 4778.23130 5289.73344 2901.128629 1513.90536 5039.378029
## 9 11209.52075 4778.23130 5289.73344 2901.128629 1513.90536 5039.378029
## 10 11209.52075 4778.23130 5289.73344 2901.128629 1513.90536 5039.378029
## 11 90.11499 46.92207 108.99533 58.068416 30.91458 40.930937
## 12 441.75934 158.36200 130.62138 83.282334 64.55693 277.745647
## 13 441.75934 158.36200 130.62138 83.282334 64.55693 277.745647
## 14 441.75934 158.36200 130.62138 83.282334 64.55693 277.745647
## 15 441.75934 158.36200 130.62138 83.282334 64.55693 277.745647
## 16 441.75934 158.36200 130.62138 83.282334 64.55693 277.745647
## 17 441.75934 158.36200 130.62138 83.282334 64.55693 277.745647
## 18 441.75934 158.36200 130.62138 83.282334 64.55693 277.745647
## RS11D RS12A RS12C RS13A RS13C RS14B RS14C
## 1 30.365093 0.00000 0.00000 52.847936 48.61116 32.08188 34.29794
## 2 7.591273 0.00000 0.00000 2.935996 5.71896 0.00000 14.07095
## 3 7.591273 0.00000 0.00000 2.935996 5.71896 0.00000 14.07095
## 4 7.591273 0.00000 0.00000 2.935996 5.71896 0.00000 14.07095
## 5 18.435949 0.00000 15.75948 21.530641 30.50112 10.93701 27.26247
## 6 600.795059 120.41737 140.08424 640.047229 997.00531 311.34011 410.69585
## 7 2731.773922 93.88473 507.80536 4233.706903 5611.25264 1408.68638 2244.31651
## 8 2731.773922 93.88473 507.80536 4233.706903 5611.25264 1408.68638 2244.31651
## 9 2731.773922 93.88473 507.80536 4233.706903 5611.25264 1408.68638 2244.31651
## 10 2731.773922 93.88473 507.80536 4233.706903 5611.25264 1408.68638 2244.31651
## 11 37.956367 56.12674 49.90501 30.338630 61.95540 37.91495 42.21285
## 12 156.163337 26.53264 35.02106 99.823880 154.41191 64.89290 89.70231
## 13 156.163337 26.53264 35.02106 99.823880 154.41191 64.89290 89.70231
## 14 156.163337 26.53264 35.02106 99.823880 154.41191 64.89290 89.70231
## 15 156.163337 26.53264 35.02106 99.823880 154.41191 64.89290 89.70231
## 16 156.163337 26.53264 35.02106 99.823880 154.41191 64.89290 89.70231
## 17 156.163337 26.53264 35.02106 99.823880 154.41191 64.89290 89.70231
## 18 156.163337 26.53264 35.02106 99.823880 154.41191 64.89290 89.70231
## RS15B RS15D RS1B RS1C RS2B RS2C
## 1 35.264808 17.090512 17.058689 29.758633 46.618913 29.996748
## 2 6.880938 0.000000 0.000000 5.759735 4.398011 9.229769
## 3 6.880938 0.000000 0.000000 5.759735 4.398011 9.229769
## 4 6.880938 0.000000 0.000000 5.759735 4.398011 9.229769
## 5 22.363049 9.202583 7.960722 0.000000 0.000000 5.768605
## 6 455.002031 118.318928 377.565660 435.819977 494.336402 485.716581
## 7 2819.464384 420.689522 2261.982225 2923.065702 3720.717045 1863.259569
## 8 2819.464384 420.689522 2261.982225 2923.065702 3720.717045 1863.259569
## 9 2819.464384 420.689522 2261.982225 2923.065702 3720.717045 1863.259569
## 10 2819.464384 420.689522 2261.982225 2923.065702 3720.717045 1863.259569
## 11 67.089146 18.405167 19.333181 62.397133 51.896526 39.226517
## 12 131.597941 60.474119 65.960266 160.312635 147.773159 122.294436
## 13 131.597941 60.474119 65.960266 160.312635 147.773159 122.294436
## 14 131.597941 60.474119 65.960266 160.312635 147.773159 122.294436
## 15 131.597941 60.474119 65.960266 160.312635 147.773159 122.294436
## 16 131.597941 60.474119 65.960266 160.312635 147.773159 122.294436
## 17 131.597941 60.474119 65.960266 160.312635 147.773159 122.294436
## 18 131.597941 60.474119 65.960266 160.312635 147.773159 122.294436
## RS3B RS3D RS6A RS6D RS7B RS7C
## 1 44.726494 27.82347 54.02640 43.647109 21.99057 40.041998
## 2 8.750836 3.83772 0.00000 12.774764 0.00000 20.020999
## 3 8.750836 3.83772 0.00000 12.774764 0.00000 20.020999
## 4 8.750836 3.83772 0.00000 12.774764 0.00000 20.020999
## 5 16.529356 20.14803 23.15417 8.516509 28.58774 7.280363
## 6 711.734643 413.51432 774.89295 465.214309 417.82076 392.229567
## 7 3845.506163 3217.92818 3720.10359 2239.841892 1795.52976 2086.734101
## 8 3845.506163 3217.92818 3720.10359 2239.841892 1795.52976 2086.734101
## 9 3845.506163 3217.92818 3720.10359 2239.841892 1795.52976 2086.734101
## 10 3845.506163 3217.92818 3720.10359 2239.841892 1795.52976 2086.734101
## 11 52.505015 85.38927 47.85195 61.744691 51.67783 61.883087
## 12 158.487359 125.68533 216.10560 137.328709 63.77264 76.443814
## 13 158.487359 125.68533 216.10560 137.328709 63.77264 76.443814
## 14 158.487359 125.68533 216.10560 137.328709 63.77264 76.443814
## 15 158.487359 125.68533 216.10560 137.328709 63.77264 76.443814
## 16 158.487359 125.68533 216.10560 137.328709 63.77264 76.443814
## 17 158.487359 125.68533 216.10560 137.328709 63.77264 76.443814
## 18 158.487359 125.68533 216.10560 137.328709 63.77264 76.443814
## RS8B RS8C RS9A RS9C accessionnumber.geneID
## 1 18.88992 14.017577 39.13713 45.883040 AJQ31790.1
## 2 0.00000 0.000000 0.00000 6.221429 aug_v2a.06327.t1
## 3 0.00000 0.000000 0.00000 6.221429 Gene:g13552
## 4 0.00000 0.000000 0.00000 6.221429 JR972076.1
## 5 0.00000 8.410546 0.00000 27.996431 Gene:g5735.t1
## 6 281.54978 309.321204 606.62546 578.592914 XP_022801463.1
## 7 1326.79210 1195.232086 3236.87057 4895.487090 ACE95141.1
## 8 1326.79210 1195.232086 3236.87057 4895.487090 EU532164.1
## 9 1326.79210 1195.232086 3236.87057 4895.487090 Gene:g29033.t1
## 10 1326.79210 1195.232086 3236.87057 4895.487090 Gene:g29034.t1
## 11 74.66017 43.921742 54.10132 48.216076 XP_022806664.1
## 12 63.86593 73.825907 233.67167 150.869658 aug_v2a.05945.t1
## 13 63.86593 73.825907 233.67167 150.869658 Gene:g2829
## 14 63.86593 73.825907 233.67167 150.869658 Gene:g2829.t1
## 15 63.86593 73.825907 233.67167 150.869658 P14_g9951
## 16 63.86593 73.825907 233.67167 150.869658 P3_g12510
## 17 63.86593 73.825907 233.67167 150.869658 P5_g11674
## 18 63.86593 73.825907 233.67167 150.869658 XP_022783415.1
## definition
## 1 solute carrier family 4 member gamma [Stylophora pistillata]
## 2 SAARP3
## 3 Acidic SOMP (Full-Length p27)
## 4 Acidic skeletal organic matrix protein (Acidic SOMP)
## 5 Annotated: Tolloid-Like
## 6 sodium bicarbonate cotransporter 3-like isoform X2
## 7 carbonic anhydrase [Stylophora pistillata]
## 8 carbonic anhydrase 2
## 9 Annotated: Carbonic Anhydrase (STPCA2-1)
## 10 Annotated: CarbonicAnhyrase
## 11 protein lingerer-like [Stylophora pistillata]
## 12 TSP-1 and VWA domain-containing
## 13 Annotated: Thrombospondin-like protein (Thrombospondin)
## 14 Annotated: Coadhesin
## 15 clone g9951 alpha collagen-like protein gene
## 16 Thrombospondin
## 17 Hemicentin
## 18 coadhesin-like isoform X3 [Stylophora pistillata]
## Ref
## 1 Zoccola et al., 2015
## 2 Takeuchi et al., 2016
## 3 Mummadisetti et al., 2021
## 4 Ramos-Silva et al., 2013
## 5 Mummadisetti et al., 2021
## 6 Zoccola et al., 2015
## 7 Moya et al., 2008
## 8 Bertucci et al., 2011
## 9 Mummadisetti et al., 2021
## 10 Mummadisetti et al., 2021
## 11 Peled et al., 2020
## 12 Takeuchi et al., 2016
## 13 Mummadisetti et al., 2021
## 14 Mummadisetti et al., 2021
## 15 Drake et al., 2013
## 16 Drake et al., 2013
## 17 Drake et al., 2013
## 18 Peled et al., 2020red_genes_data_DE_Origin <- red_genes_data %>% filter(Origin < 0.05) %>% dplyr::select(c(Gene, Origin, definition, Ref))
red_genes_data_DE_Origin## Gene Origin
## 1 Pocillopora_acuta_HIv2___TS.g12304.t1 0.01164322
## 2 Pocillopora_acuta_HIv2___TS.g12304.t1 0.01164322
## 3 Pocillopora_acuta_HIv2___TS.g12304.t1 0.01164322
## 4 Pocillopora_acuta_HIv2___TS.g12304.t1 0.01164322
## 5 Pocillopora_acuta_HIv2___RNAseq.g7908.t1 0.00158533
## definition Ref
## 1 carbonic anhydrase [Stylophora pistillata] Moya et al., 2008
## 2 carbonic anhydrase 2 Bertucci et al., 2011
## 3 Annotated: Carbonic Anhydrase (STPCA2-1) Mummadisetti et al., 2021
## 4 Annotated: CarbonicAnhyrase Mummadisetti et al., 2021
## 5 protein lingerer-like [Stylophora pistillata] Peled et al., 2020unique(red_genes_data_DE_Origin$Gene)## [1] "Pocillopora_acuta_HIv2___TS.g12304.t1"
## [2] "Pocillopora_acuta_HIv2___RNAseq.g7908.t1"Carbonic Anhydrase (STPCA2-1)
biomin_all_filtered_long_g12304 <- biomin_all_filtered_long %>%
filter(Gene == "Pocillopora_acuta_HIv2___TS.g12304.t1")
biomin_all_filtered_long_g12304 ## # A tibble: 184 × 36
## Gene Dispersion AIC logLik meanExp X.Intercept. TreatmentVariable
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Pocillopora_a… 0.457 841. -415. 11.6 7.88 0.0448
## 2 Pocillopora_a… 0.457 841. -415. 11.6 7.88 0.0448
## 3 Pocillopora_a… 0.457 841. -415. 11.6 7.88 0.0448
## 4 Pocillopora_a… 0.457 841. -415. 11.6 7.88 0.0448
## 5 Pocillopora_a… 0.457 841. -415. 11.6 7.88 0.0448
## 6 Pocillopora_a… 0.457 841. -415. 11.6 7.88 0.0448
## 7 Pocillopora_a… 0.457 841. -415. 11.6 7.88 0.0448
## 8 Pocillopora_a… 0.457 841. -415. 11.6 7.88 0.0448
## 9 Pocillopora_a… 0.457 841. -415. 11.6 7.88 0.0448
## 10 Pocillopora_a… 0.457 841. -415. 11.6 7.88 0.0448
## # ℹ 174 more rows
## # ℹ 29 more variables: OriginFlat <dbl>, TreatmentVariable.OriginFlat <dbl>,
## # se_.Intercept. <dbl>, se_TreatmentVariable <dbl>, se_OriginFlat <dbl>,
## # se_TreatmentVariable.OriginFlat <dbl>, Chisq_Treatment <dbl>,
## # Chisq_Origin <dbl>, Chisq_Treatment.Origin <dbl>, Df_Treatment <int>,
## # Df_Origin <int>, Df_Treatment.Origin <int>, P_Treatment <dbl>,
## # P_Origin <dbl>, P_Treatment.Origin <dbl>, Treatment <fct>, …library(nlme)##
## Attaching package: 'nlme'## The following object is masked from 'package:IRanges':
##
## collapse## The following object is masked from 'package:dplyr':
##
## collapselibrary(emmeans)
g12304.lme <- lme(Expression~Origin*Treatment2, random = ~1|Colony, data=biomin_all_filtered_long_g12304 , na.action=na.exclude)
car::Anova(g12304.lme, type=3)## Analysis of Deviance Table (Type III tests)
##
## Response: Expression
## Chisq Df Pr(>Chisq)
## (Intercept) 196.964 1 < 2.2e-16 ***
## Origin 55.461 1 9.531e-14 ***
## Treatment2 14.008 1 0.0001820 ***
## Origin:Treatment2 11.399 1 0.0007348 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1tukey3<- emmeans(g12304.lme, list(pairwise ~ Origin), adjust = "tukey")## NOTE: Results may be misleading due to involvement in interactionstukey3## $`emmeans of Origin`
## Origin emmean SE df lower.CL upper.CL
## RF 5238 384 19 4434 6042
## RS 2770 375 19 1985 3555
##
## Results are averaged over the levels of: Treatment2
## Degrees-of-freedom method: containment
## Confidence level used: 0.95
##
## $`pairwise differences of Origin`
## 1 estimate SE df t.ratio p.value
## RF - RS 2468 371 161 6.657 <.0001
##
## Results are averaged over the levels of: Treatment2
## Degrees-of-freedom method: containmentlibrary(Rmisc)
g12304_sum<-summarySE(biomin_all_filtered_long_g12304 , measurevar='Expression', groupvars=c('Origin', 'Treatment2'), na.rm=TRUE, conf.interval = 0.95)
g12304_sum## Origin Treatment2 N Expression sd se ci
## 1 RF Stable 44 5960.191 2423.309 365.3276 736.7533
## 2 RF Variable 44 4766.484 2433.696 366.8934 739.9111
## 3 RS Stable 48 2494.782 1512.223 218.2706 439.1038
## 4 RS Variable 48 2791.885 1392.250 200.9540 404.2674Figure
pd<- position_dodge(0.2)
g12304_fig<-ggplot(data=g12304_sum, aes(y=Expression, x=Treatment2, color=Origin),group = interaction(Origin))+
geom_point(data=biomin_all_filtered_long_g12304,aes(y=Expression, x=Treatment2, color=Origin), alpha=0.4, position = pd)+
geom_line(aes(group = interaction(Origin), stat="identity"),position=position_dodge(0.2))+
geom_point(size=3, stat="identity", position = pd)+
geom_errorbar(aes(ymin=Expression-se, ymax=Expression+se), stat="identity",width=0.2, position = pd)+
scale_fill_manual("Origin", values=c("RS"='#26519B', "RF"= "#FE180C"))+
scale_color_manual("Origin", values=c("RS"='#26519B', "RF"= "#FE180C"))+
scale_y_continuous(expression(Carbonic~Anhydrase~(STPCA2~1)))+
ggtitle(~MERed)+
theme_classic()+
theme(axis.text.x=element_text(vjust=0.5,size=12),#angling the labels on the x-axis
plot.title = element_text(margin = margin(t = 10, b = 10), hjust=0.5),#telling it where to position our plot title
panel.background= element_rect(fill=NA, color='black'),#this is making the black box around the graph
#strip.background = element_blank(),
#strip.text = element_blank(),
legend.title = element_text(vjust=0.5,size=12),
legend.position="none",
axis.text.y = element_text(vjust=0.5, size=12), #making the axis text larger
axis.title.x = element_blank(),#making the axis title larger
axis.title.y = element_text(size=12))#making the axis title larger ## Warning in geom_line(aes(group = interaction(Origin), stat = "identity"), :
## Ignoring unknown aesthetics: statg12304_fig
Protein lingerer-like
biomin_all_filtered_long_g7908 <- biomin_all_filtered_long %>%
filter(Gene == "Pocillopora_acuta_HIv2___RNAseq.g7908.t1")
biomin_all_filtered_long_g7908 ## # A tibble: 46 × 36
## Gene Dispersion AIC logLik meanExp X.Intercept. TreatmentVariable
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Pocillopora_a… 0.138 426. -207. 5.95 3.98 -0.0743
## 2 Pocillopora_a… 0.138 426. -207. 5.95 3.98 -0.0743
## 3 Pocillopora_a… 0.138 426. -207. 5.95 3.98 -0.0743
## 4 Pocillopora_a… 0.138 426. -207. 5.95 3.98 -0.0743
## 5 Pocillopora_a… 0.138 426. -207. 5.95 3.98 -0.0743
## 6 Pocillopora_a… 0.138 426. -207. 5.95 3.98 -0.0743
## 7 Pocillopora_a… 0.138 426. -207. 5.95 3.98 -0.0743
## 8 Pocillopora_a… 0.138 426. -207. 5.95 3.98 -0.0743
## 9 Pocillopora_a… 0.138 426. -207. 5.95 3.98 -0.0743
## 10 Pocillopora_a… 0.138 426. -207. 5.95 3.98 -0.0743
## # ℹ 36 more rows
## # ℹ 29 more variables: OriginFlat <dbl>, TreatmentVariable.OriginFlat <dbl>,
## # se_.Intercept. <dbl>, se_TreatmentVariable <dbl>, se_OriginFlat <dbl>,
## # se_TreatmentVariable.OriginFlat <dbl>, Chisq_Treatment <dbl>,
## # Chisq_Origin <dbl>, Chisq_Treatment.Origin <dbl>, Df_Treatment <int>,
## # Df_Origin <int>, Df_Treatment.Origin <int>, P_Treatment <dbl>,
## # P_Origin <dbl>, P_Treatment.Origin <dbl>, Treatment <fct>, …library(nlme)
library(emmeans)
g7908.lme <- lme(Expression~Origin*Treatment2, random = ~1|Colony, data=biomin_all_filtered_long_g7908 , na.action=na.exclude)
car::Anova(g7908.lme, type=3)## Analysis of Deviance Table (Type III tests)
##
## Response: Expression
## Chisq Df Pr(>Chisq)
## (Intercept) 207.4212 1 < 2.2e-16 ***
## Origin 15.8613 1 6.816e-05 ***
## Treatment2 0.1605 1 0.6887
## Origin:Treatment2 0.0058 1 0.9392
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1tukey3<- emmeans(g7908.lme, list(pairwise ~ Origin), adjust = "tukey")## NOTE: Results may be misleading due to involvement in interactionstukey3## $`emmeans of Origin`
## Origin emmean SE df lower.CL upper.CL
## RF 81.2 4.07 19 72.7 89.7
## RS 49.9 3.89 19 41.8 58.0
##
## Results are averaged over the levels of: Treatment2
## Degrees-of-freedom method: containment
## Confidence level used: 0.95
##
## $`pairwise differences of Origin`
## 1 estimate SE df t.ratio p.value
## RF - RS 31.3 5.63 23 5.556 <.0001
##
## Results are averaged over the levels of: Treatment2
## Degrees-of-freedom method: containmentlibrary(Rmisc)
g7908_sum<-summarySE(biomin_all_filtered_long_g7908 , measurevar='Expression', groupvars=c('Origin', 'Treatment2'), na.rm=TRUE, conf.interval = 0.95)
g7908_sum## Origin Treatment2 N Expression sd se ci
## 1 RF Stable 11 82.80022 22.33628 6.734642 15.005717
## 2 RF Variable 11 79.54283 21.44921 6.467180 14.409774
## 3 RS Stable 12 51.10111 16.98594 4.903419 10.792352
## 4 RS Variable 12 48.70220 15.09645 4.357969 9.591824Figure
pd<- position_dodge(0.2)
g7908_fig<-ggplot(data=g7908_sum, aes(y=Expression, x=Treatment2, color=Origin),group = interaction(Origin))+
geom_point(data=biomin_all_filtered_long_g7908,aes(y=Expression, x=Treatment2, color=Origin), alpha=0.4, position = pd)+
geom_line(aes(group = interaction(Origin), stat="identity"),position=position_dodge(0.2))+
geom_point(size=3, stat="identity", position = pd)+
geom_errorbar(aes(ymin=Expression-se, ymax=Expression+se), stat="identity",width=0.2, position = pd)+
scale_fill_manual("Origin", values=c("RS"='#26519B', "RF"= "#FE180C"))+
scale_color_manual("Origin", values=c("RS"='#26519B', "RF"= "#FE180C"))+
scale_y_continuous(expression(Protein~lingerer-like))+
ggtitle(~MERed)+
theme_classic()+
theme(axis.text.x=element_text(vjust=0.5,size=12),#angling the labels on the x-axis
plot.title = element_text(margin = margin(t = 10, b = 10), hjust=0.5),#telling it where to position our plot title
panel.background= element_rect(fill=NA, color='black'),#this is making the black box around the graph
#strip.background = element_blank(),
#strip.text = element_blank(),
legend.title = element_text(vjust=0.5,size=12),
legend.position="none",
axis.text.y = element_text(vjust=0.5, size=12), #making the axis text larger
axis.title.x = element_blank(),#making the axis title larger
axis.title.y = element_text(size=12))#making the axis title larger ## Warning in geom_line(aes(group = interaction(Origin), stat = "identity"), :
## Ignoring unknown aesthetics: statg7908_fig
MEBrown signficant genes
brown_genes <- biomin_mod %>% filter(moduleColor == "brown") %>% pull(Pocillopora_acuta_best_hit) %>% unique()
length(brown_genes)## [1] 10brown_genes_data <- biomin_all %>% filter(Gene %in% brown_genes)
brown_genes_data## Gene Dispersion AIC logLik
## 1 Pocillopora_acuta_HIv2___RNAseq.g25351.t1 0.13435721 815.9720 -401.9860
## 2 Pocillopora_acuta_HIv2___TS.g13222.t1b 0.05064785 621.2754 -304.6377
## 3 Pocillopora_acuta_HIv2___TS.g13222.t1b 0.05064785 621.2754 -304.6377
## 4 Pocillopora_acuta_HIv2___TS.g13222.t1b 0.05064785 621.2754 -304.6377
## 5 Pocillopora_acuta_HIv2___TS.g13222.t1b 0.05064785 621.2754 -304.6377
## 6 Pocillopora_acuta_HIv2___RNAseq.g21232.t1 0.06053790 353.3934 -170.6967
## 7 Pocillopora_acuta_HIv2___RNAseq.g5013.t1 0.47107150 304.1954 -146.0977
## 8 Pocillopora_acuta_HIv2___RNAseq.g10093.t2 0.08759845 385.1631 -186.5816
## 9 Pocillopora_acuta_HIv2___RNAseq.g13824.t1 2.35237053 452.5561 -220.2780
## 10 Pocillopora_acuta_HIv2___TS.g425.t1 3.08083166 263.7848 -125.8924
## 11 Pocillopora_acuta_HIv2___RNAseq.g22388.t1 0.02224702 453.8740 -220.9370
## 12 Pocillopora_acuta_HIv2___RNAseq.g22388.t1 0.02224702 453.8740 -220.9370
## 13 Pocillopora_acuta_HIv2___TS.g5338.t1 0.30589381 351.9165 -169.9583
## 14 Pocillopora_acuta_HIv2___RNAseq.g22261.t1 0.49716428 326.8663 -157.4332
## 15 Pocillopora_acuta_HIv2___RNAseq.g22261.t1 0.49716428 326.8663 -157.4332
## 16 Pocillopora_acuta_HIv2___RNAseq.g22261.t1 0.49716428 326.8663 -157.4332
## 17 Pocillopora_acuta_HIv2___RNAseq.g22261.t1 0.49716428 326.8663 -157.4332
## meanExp X.Intercept. TreatmentVariable OriginFlat
## 1 12.237463 8.288576 -0.017091983 0.4044018
## 2 9.498170 6.531267 -0.123332597 0.2246403
## 3 9.498170 6.531267 -0.123332597 0.2246403
## 4 9.498170 6.531267 -0.123332597 0.2246403
## 5 9.498170 6.531267 -0.123332597 0.2246403
## 6 5.221702 3.540412 -0.060041477 0.2071984
## 7 3.161833 1.828096 0.206488581 0.7336259
## 8 5.739306 3.881564 -0.141912018 0.2630129
## 9 4.663709 2.363210 -0.064800069 2.6486369
## 10 1.964418 1.517871 -1.214684460 0.6328338
## 11 7.225058 4.954065 0.004692091 0.1021809
## 12 7.225058 4.954065 0.004692091 0.1021809
## 13 4.057441 2.839322 0.004254032 0.1654600
## 14 3.399138 2.360264 -0.048347784 0.3052264
## 15 3.399138 2.360264 -0.048347784 0.3052264
## 16 3.399138 2.360264 -0.048347784 0.3052264
## 17 3.399138 2.360264 -0.048347784 0.3052264
## TreatmentVariable.OriginFlat se_.Intercept. se_TreatmentVariable
## 1 0.22170889 0.10590424 0.14977519
## 2 0.12611496 0.07177003 0.09372265
## 3 0.12611496 0.07177003 0.09372265
## 4 0.12611496 0.07177003 0.09372265
## 5 0.12611496 0.07177003 0.09372265
## 6 0.07052439 0.08714127 0.12284578
## 7 -0.09755882 0.25241984 0.32039420
## 8 0.21951326 0.09496357 0.13527058
## 9 -0.06976957 0.45152756 0.63897743
## 10 1.43782801 0.52440626 0.77024650
## 11 0.06248481 0.04941464 0.06986315
## 12 0.06248481 0.04941464 0.06986315
## 13 0.03980596 0.17425056 0.24638553
## 14 0.09166469 0.22202875 0.31461588
## 15 0.09166469 0.22202875 0.31461588
## 16 0.09166469 0.22202875 0.31461588
## 17 0.09166469 0.22202875 0.31461588
## se_OriginFlat se_TreatmentVariable.OriginFlat Chisq_Treatment Chisq_Origin
## 1 0.15311285 0.2165347 0.6766804507 22.648367
## 2 0.10247292 0.1351560 0.8617790376 14.103484
## 3 0.10247292 0.1351560 0.8617790376 14.103484
## 4 0.10247292 0.1351560 0.8617790376 14.103484
## 5 0.10247292 0.1351560 0.8617790376 14.103484
## 6 0.12365034 0.1742061 0.0822376051 7.566016
## 7 0.33365787 0.4499480 0.4877012204 7.766492
## 8 0.13572952 0.1923522 0.1202652502 14.986002
## 9 0.64679103 0.9150541 0.0466828440 32.638688
## 10 0.75210323 1.0828096 0.8096894626 6.010551
## 11 0.07101623 0.1002309 0.4895241958 7.101397
## 12 0.07101623 0.1002309 0.4895241958 7.101397
## 13 0.25035235 0.3537481 0.0177651612 1.098700
## 14 0.31751820 0.4491738 0.0002261147 2.443000
## 15 0.31751820 0.4491738 0.0002261147 2.443000
## 16 0.31751820 0.4491738 0.0002261147 2.443000
## 17 0.31751820 0.4491738 0.0002261147 2.443000
## Chisq_Treatment.Origin Df_Treatment Df_Origin Df_Treatment.Origin
## 1 1.04836235 1 1 1
## 2 0.87068820 1 1 1
## 3 0.87068820 1 1 1
## 4 0.87068820 1 1 1
## 5 0.87068820 1 1 1
## 6 0.16388974 1 1 1
## 7 0.04701196 1 1 1
## 8 1.30234851 1 1 1
## 9 0.00581351 1 1 1
## 10 1.76323289 1 1 1
## 11 0.38863841 1 1 1
## 12 0.38863841 1 1 1
## 13 0.01266216 1 1 1
## 14 0.04164618 1 1 1
## 15 0.04164618 1 1 1
## 16 0.04164618 1 1 1
## 17 0.04164618 1 1 1
## P_Treatment P_Origin P_Treatment.Origin Treatment Origin
## 1 0.4107321 1.945255e-06 0.3058846 1 1.683701e-04
## 2 0.3532413 1.730230e-04 0.3507649 1 4.271694e-03
## 3 0.3532413 1.730230e-04 0.3507649 1 4.271694e-03
## 4 0.3532413 1.730230e-04 0.3507649 1 4.271694e-03
## 5 0.3532413 1.730230e-04 0.3507649 1 4.271694e-03
## 6 0.7742877 5.947912e-03 0.6856003 1 4.567166e-02
## 7 0.4849545 5.322430e-03 0.8283467 1 4.281486e-02
## 8 0.7287470 1.083116e-04 0.2537847 1 3.022284e-03
## 9 0.8289393 1.109835e-08 0.9392231 1 2.246821e-06
## 10 0.3682121 1.422058e-02 0.1842218 1 7.714976e-02
## 11 0.4841396 7.702390e-03 0.5330160 1 5.303815e-02
## 12 0.4841396 7.702390e-03 0.5330160 1 5.303815e-02
## 13 0.8939672 2.945516e-01 0.9104061 1 4.212215e-01
## 14 0.9880026 1.180503e-01 0.8382957 1 2.665372e-01
## 15 0.9880026 1.180503e-01 0.8382957 1 2.665372e-01
## 16 0.9880026 1.180503e-01 0.8382957 1 2.665372e-01
## 17 0.9880026 1.180503e-01 0.8382957 1 2.665372e-01
## Treatment.Origin Stable_OriginFC Variable_OriginFC maxGroupFC
## 1 0.9994279 -0.5833078 -0.9031151 Variable
## 2 0.9994279 -0.3236650 -0.5053302 Variable
## 3 0.9994279 -0.3236650 -0.5053302 Variable
## 4 0.9994279 -0.3236650 -0.5053302 Variable
## 5 0.9994279 -0.3236650 -0.5053302 Variable
## 6 0.9994279 -0.2912946 -0.3901793 Variable
## 7 0.9994279 -0.9506277 -0.8369316 Stable
## 8 0.9994279 -0.3726912 -0.6832631 Variable
## 9 0.9994279 -3.7009819 -3.5934153 Stable
## 10 0.9994279 -0.7859341 -2.3179765 Variable
## 11 0.9994279 -0.1464338 -0.2360343 Variable
## 12 0.9994279 -0.1464338 -0.2360343 Variable
## 13 0.9994279 -0.2265047 -0.2813321 Variable
## 14 0.9994279 -0.4072391 -0.5293491 Variable
## 15 0.9994279 -0.4072391 -0.5293491 Variable
## 16 0.9994279 -0.4072391 -0.5293491 Variable
## 17 0.9994279 -0.4072391 -0.5293491 Variable
## col RF13B RF13D RF14B RF14C RF15B
## 1 q_Origin < 0.05 5832.107578 4927.06132 13021.916529 6113.266917 5054.03248
## 2 q_Origin < 0.05 768.436933 898.68343 641.209339 1029.726405 929.73223
## 3 q_Origin < 0.05 768.436933 898.68343 641.209339 1029.726405 929.73223
## 4 q_Origin < 0.05 768.436933 898.68343 641.209339 1029.726405 929.73223
## 5 q_Origin < 0.05 768.436933 898.68343 641.209339 1029.726405 929.73223
## 6 q_Origin < 0.05 22.601086 27.47067 40.487521 27.743924 57.47802
## 7 q_Origin < 0.05 6.163933 13.73534 13.181983 6.402444 23.19289
## 8 q_Origin < 0.05 60.612004 67.69559 80.975041 80.030550 49.41093
## 9 q_Origin < 0.05 142.797772 160.89965 159.125372 265.701425 175.45923
## 10 Not Significant 3.081966 0.00000 10.357273 0.000000 16.13418
## 11 Not Significant 136.633840 135.39117 153.475950 166.463543 165.37536
## 12 Not Significant 136.633840 135.39117 153.475950 166.463543 165.37536
## 13 Not Significant 15.409832 13.73534 8.474132 37.347590 22.18450
## 14 Not Significant 0.000000 0.00000 16.948264 13.871962 13.10902
## 15 Not Significant 0.000000 0.00000 16.948264 13.871962 13.10902
## 16 Not Significant 0.000000 0.00000 16.948264 13.871962 13.10902
## 17 Not Significant 0.000000 0.00000 16.948264 13.871962 13.10902
## RF15D RF17B RF17D RF18B RF18D RF19B
## 1 7484.794135 8321.89749 7240.96832 5891.93841 5089.20850 8016.94162
## 2 729.228491 662.64670 829.34143 873.57556 1006.69338 874.58961
## 3 729.228491 662.64670 829.34143 873.57556 1006.69338 874.58961
## 4 729.228491 662.64670 829.34143 873.57556 1006.69338 874.58961
## 5 729.228491 662.64670 829.34143 873.57556 1006.69338 874.58961
## 6 52.341218 46.12651 47.89155 59.88220 56.85644 34.26709
## 7 7.984254 18.00059 24.52982 25.83154 13.37799 17.91234
## 8 55.889775 90.00295 66.58093 115.06775 81.38274 74.76456
## 9 205.816314 261.00855 238.28965 286.49521 376.81325 59.96741
## 10 0.000000 14.62548 24.52982 12.91577 11.14832 10.90317
## 11 160.572210 146.25479 216.09601 184.34323 160.53582 165.10507
## 12 160.572210 146.25479 216.09601 184.34323 160.53582 165.10507
## 13 14.194229 27.00088 35.04260 24.65737 24.52631 14.79715
## 14 14.194229 21.37570 22.19364 25.83154 28.98563 18.69114
## 15 14.194229 21.37570 22.19364 25.83154 28.98563 18.69114
## 16 14.194229 21.37570 22.19364 25.83154 28.98563 18.69114
## 17 14.194229 21.37570 22.19364 25.83154 28.98563 18.69114
## RF19C RF20B RF20C RF22B RF22C RF23A
## 1 8866.41520 4982.125947 8081.53904 7060.316045 3512.969287 6530.045150
## 2 809.07623 788.720960 1191.92433 840.090197 774.146857 991.744488
## 3 809.07623 788.720960 1191.92433 840.090197 774.146857 991.744488
## 4 809.07623 788.720960 1191.92433 840.090197 774.146857 991.744488
## 5 809.07623 788.720960 1191.92433 840.090197 774.146857 991.744488
## 6 43.79615 37.007435 35.59345 37.904622 50.531779 27.412483
## 7 26.50820 9.251859 21.90366 13.150583 2.021271 8.811155
## 8 89.89736 47.415776 69.36160 58.790842 53.563686 46.992829
## 9 165.96436 188.506622 86.70200 37.904622 91.967838 57.762019
## 10 0.00000 10.408341 14.60244 15.471274 14.148898 0.000000
## 11 150.98146 232.452952 135.07259 152.392052 168.776142 140.978486
## 12 150.98146 232.452952 135.07259 152.392052 168.776142 140.978486
## 13 20.74554 45.102811 11.86448 13.924147 18.191440 17.622311
## 14 13.83036 39.320400 20.07836 5.414946 11.116991 11.748207
## 15 13.83036 39.320400 20.07836 5.414946 11.116991 11.748207
## 16 13.83036 39.320400 20.07836 5.414946 11.116991 11.748207
## 17 13.83036 39.320400 20.07836 5.414946 11.116991 11.748207
## RF23C RF24B RF24D RF25A RF25C RS11B
## 1 4346.08909 4820.265658 5124.510370 6260.692121 4758.11817 4408.846684
## 2 775.77250 963.857623 592.553951 705.225632 818.32722 602.269507
## 3 775.77250 963.857623 592.553951 705.225632 818.32722 602.269507
## 4 775.77250 963.857623 592.553951 705.225632 818.32722 602.269507
## 5 775.77250 963.857623 592.553951 705.225632 818.32722 602.269507
## 6 39.18043 38.124186 50.172451 50.427835 33.64234 33.134568
## 7 11.75413 11.730519 11.245549 7.640581 7.27402 1.949092
## 8 41.13945 56.697507 38.061860 45.079428 58.19216 36.058207
## 9 86.19694 79.181001 21.626057 43.551312 20.00355 13.643646
## 10 20.56973 8.797889 3.460169 9.168697 7.27402 0.000000
## 11 171.41438 146.631484 151.382396 149.755389 115.47506 180.291034
## 12 171.41438 146.631484 151.382396 149.755389 115.47506 180.291034
## 13 24.48777 16.618235 11.245549 22.921743 16.36654 13.643646
## 14 13.71315 7.820346 14.705718 6.112465 7.27402 14.618192
## 15 13.71315 7.820346 14.705718 6.112465 7.27402 14.618192
## 16 13.71315 7.820346 14.705718 6.112465 7.27402 14.618192
## 17 13.71315 7.820346 14.705718 6.112465 7.27402 14.618192
## RS11D RS12A RS12C RS13A RS13C RS14B
## 1 3804.312393 2168.533115 3681.588844 2701.116748 2417.21364 3177.564822
## 2 638.751425 600.045869 607.615376 652.769881 770.15324 740.070742
## 3 638.751425 600.045869 607.615376 652.769881 770.15324 740.070742
## 4 638.751425 600.045869 607.615376 652.769881 770.15324 740.070742
## 5 638.751425 600.045869 607.615376 652.769881 770.15324 740.070742
## 6 34.702964 14.286806 26.265794 29.359965 24.78216 29.894483
## 7 5.422338 3.061459 7.879738 19.573310 19.06320 6.562204
## 8 48.801043 24.491668 33.270006 36.210623 33.36060 32.081884
## 9 0.000000 0.000000 0.000000 9.786655 14.29740 0.000000
## 10 6.506806 0.000000 6.128685 8.807989 0.00000 0.000000
## 11 125.798243 156.134384 138.333184 132.119841 145.83347 99.162188
## 12 125.798243 156.134384 138.333184 132.119841 145.83347 99.162188
## 13 10.844676 0.000000 0.000000 27.402634 34.31376 11.666140
## 14 7.591273 2.040972 3.502106 0.000000 18.11004 8.749605
## 15 7.591273 2.040972 3.502106 0.000000 18.11004 8.749605
## 16 7.591273 2.040972 3.502106 0.000000 18.11004 8.749605
## 17 7.591273 2.040972 3.502106 0.000000 18.11004 8.749605
## RS14C RS15B RS15D RS1B RS1C RS2B
## 1 3123.75088 2933.859980 3892.692731 3421.973109 3717.909183 6782.61209
## 2 784.45546 556.495868 402.284355 536.780096 536.615346 469.70754
## 3 784.45546 556.495868 402.284355 536.780096 536.615346 469.70754
## 4 784.45546 556.495868 402.284355 536.780096 536.615346 469.70754
## 5 784.45546 556.495868 402.284355 536.780096 536.615346 469.70754
## 6 32.53907 29.243987 35.495678 30.705641 34.558412 32.54528
## 7 0.00000 2.580352 3.943964 2.274492 1.919912 12.31443
## 8 38.69511 36.985042 57.844809 44.352593 40.318148 66.84976
## 9 0.00000 0.000000 5.258619 0.000000 23.998897 19.35125
## 10 0.00000 0.000000 0.000000 0.000000 14.399338 0.00000
## 11 148.62441 142.779465 136.724095 142.155746 113.274796 157.44878
## 12 148.62441 142.779465 136.724095 142.155746 113.274796 157.44878
## 13 0.00000 19.782697 27.607750 9.097968 21.119030 14.07363
## 14 15.82982 6.020821 9.202583 20.470427 4.799779 15.83284
## 15 15.82982 6.020821 9.202583 20.470427 4.799779 15.83284
## 16 15.82982 6.020821 9.202583 20.470427 4.799779 15.83284
## 17 15.82982 6.020821 9.202583 20.470427 4.799779 15.83284
## RS2C RS3B RS3D RS6A RS6D RS7B
## 1 5170.977950 5164.937734 4277.13888 5058.414714 4024.050548 2932.442053
## 2 625.316834 653.395737 607.31918 523.284281 498.215782 657.517941
## 3 625.316834 653.395737 607.31918 523.284281 498.215782 657.517941
## 4 625.316834 653.395737 607.31918 523.284281 498.215782 657.517941
## 5 625.316834 653.395737 607.31918 523.284281 498.215782 657.517941
## 6 39.226517 20.418617 44.13378 38.590286 38.324291 58.275001
## 7 3.461163 1.944630 1.91886 9.261669 9.581073 4.398113
## 8 65.762102 43.754179 48.93093 49.395566 37.259727 41.782076
## 9 20.766980 18.473987 6.71601 13.892503 24.484964 17.592453
## 10 0.000000 0.000000 11.51316 0.000000 0.000000 6.597170
## 11 104.988620 148.764208 153.50880 114.227247 153.297164 136.341513
## 12 104.988620 148.764208 153.50880 114.227247 153.297164 136.341513
## 13 14.998374 24.307877 21.10746 30.872229 25.549527 18.691982
## 14 27.689306 6.806206 12.47259 9.261669 8.516509 12.094812
## 15 27.689306 6.806206 12.47259 9.261669 8.516509 12.094812
## 16 27.689306 6.806206 12.47259 9.261669 8.516509 12.094812
## 17 27.689306 6.806206 12.47259 9.261669 8.516509 12.094812
## RS7C RS8B RS8C RS9A RS9C
## 1 5029.820923 3791.47710 4200.600647 4097.88734 3391.456585
## 2 796.289725 577.49189 694.337326 627.34511 960.433131
## 3 796.289725 577.49189 694.337326 627.34511 960.433131
## 4 796.289725 577.49189 694.337326 627.34511 960.433131
## 5 796.289725 577.49189 694.337326 627.34511 960.433131
## 6 44.592225 33.28224 24.297134 42.59040 27.218753
## 7 3.640182 18.88992 11.214062 10.35983 6.221429
## 8 50.052497 43.17696 52.332288 48.34586 62.214292
## 9 0.000000 11.69376 12.148567 18.41747 17.108930
## 10 2.730136 0.00000 0.000000 0.00000 8.554465
## 11 122.856129 139.42561 114.009628 142.73540 179.643768
## 12 122.856129 139.42561 114.009628 142.73540 179.643768
## 13 20.020999 19.78944 15.886588 19.56856 13.998216
## 14 0.000000 11.69376 9.345051 13.81310 9.332144
## 15 0.000000 11.69376 9.345051 13.81310 9.332144
## 16 0.000000 11.69376 9.345051 13.81310 9.332144
## 17 0.000000 11.69376 9.345051 13.81310 9.332144
## accessionnumber.geneID
## 1 XP_022794351.1
## 2 Gene:g15294.t1
## 3 P24_g15888
## 4 P26_g1441
## 5 XP_022779720.1
## 6 XP_022783044.1
## 7 JR986059.1
## 8 XP_022804785.1
## 9 Gene:g27814
## 10 P22_g19762
## 11 Gene:g24177
## 12 P23_g1057
## 13 XP_022783952.1
## 14 aug_v2a.01440.t1(aug_v2a.01441.t1)
## 15 JR991407.1
## 16 P15_g1532
## 17 XP_022780690.1
## definition
## 1 mammalian ependymin-related protein 1-like [Stylophora pistillata]
## 2 Annotated: Vitellogenin
## 3 clone g15888 vitellogenin-like protein gene
## 4 clone g1441 vitellogenin-like protein gene
## 5 vitellogenin-like [Stylophora pistillata]
## 6 uncharacterized protein LOC111323869 [Stylophora pistillata]
## 7 Cephalotoxin-like protein
## 8 thioredoxin reductase 1, cytoplasmic-like [Stylophora pistillata]
## 9 Annotated: carbonic anhydrase (STPCA2-2)
## 10 Stylophora pistillata clone g19762 hypothetical protein gene
## 11 Annotated: Protocadherin (PC5)
## 12 Protocadherin
## 13 collagenase 3-like [Stylophora pistillata]
## 14 Adi-SAARP2
## 15 Skeletal acidic Asp-rich Protein 2 (SAARP2)
## 16 CARP9
## 17 skeletal aspartic acid-rich protein 2-like (CARP5)
## Ref
## 1 Peled et al., 2020
## 2 Mummadisetti et al., 2021
## 3 Drake et al., 2013
## 4 Drake et al., 2013
## 5 Peled et al., 2020
## 6 Peled et al., 2020
## 7 Ramos-Silva et al., 2013
## 8 Peled et al., 2020
## 9 Mummadisetti et al., 2021
## 10 Drake et al., 2013
## 11 Mummadisetti et al., 2021
## 12 Drake et al., 2013
## 13 Peled et al., 2020
## 14 Takeuchi et al., 2016
## 15 Ramos-Silva et al., 2013
## 16 Drake et al., 2013
## 17 Peled et al., 2020brown_genes_data_DE_Origin <- brown_genes_data %>% filter(Origin < 0.05) %>% dplyr::select(c(Gene, Origin, definition, Ref))
brown_genes_data_DE_Origin## Gene Origin
## 1 Pocillopora_acuta_HIv2___RNAseq.g25351.t1 1.683701e-04
## 2 Pocillopora_acuta_HIv2___TS.g13222.t1b 4.271694e-03
## 3 Pocillopora_acuta_HIv2___TS.g13222.t1b 4.271694e-03
## 4 Pocillopora_acuta_HIv2___TS.g13222.t1b 4.271694e-03
## 5 Pocillopora_acuta_HIv2___TS.g13222.t1b 4.271694e-03
## 6 Pocillopora_acuta_HIv2___RNAseq.g21232.t1 4.567166e-02
## 7 Pocillopora_acuta_HIv2___RNAseq.g5013.t1 4.281486e-02
## 8 Pocillopora_acuta_HIv2___RNAseq.g10093.t2 3.022284e-03
## 9 Pocillopora_acuta_HIv2___RNAseq.g13824.t1 2.246821e-06
## definition
## 1 mammalian ependymin-related protein 1-like [Stylophora pistillata]
## 2 Annotated: Vitellogenin
## 3 clone g15888 vitellogenin-like protein gene
## 4 clone g1441 vitellogenin-like protein gene
## 5 vitellogenin-like [Stylophora pistillata]
## 6 uncharacterized protein LOC111323869 [Stylophora pistillata]
## 7 Cephalotoxin-like protein
## 8 thioredoxin reductase 1, cytoplasmic-like [Stylophora pistillata]
## 9 Annotated: carbonic anhydrase (STPCA2-2)
## Ref
## 1 Peled et al., 2020
## 2 Mummadisetti et al., 2021
## 3 Drake et al., 2013
## 4 Drake et al., 2013
## 5 Peled et al., 2020
## 6 Peled et al., 2020
## 7 Ramos-Silva et al., 2013
## 8 Peled et al., 2020
## 9 Mummadisetti et al., 2021unique(brown_genes_data_DE_Origin$Gene)## [1] "Pocillopora_acuta_HIv2___RNAseq.g25351.t1"
## [2] "Pocillopora_acuta_HIv2___TS.g13222.t1b"
## [3] "Pocillopora_acuta_HIv2___RNAseq.g21232.t1"
## [4] "Pocillopora_acuta_HIv2___RNAseq.g5013.t1"
## [5] "Pocillopora_acuta_HIv2___RNAseq.g10093.t2"
## [6] "Pocillopora_acuta_HIv2___RNAseq.g13824.t1"mammalian ependymin-related protein 1-like
biomin_all_filtered_long_g25351 <- biomin_all_filtered_long %>%
filter(Gene == "Pocillopora_acuta_HIv2___RNAseq.g25351.t1")
biomin_all_filtered_long_g25351 ## # A tibble: 46 × 36
## Gene Dispersion AIC logLik meanExp X.Intercept. TreatmentVariable
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Pocillopora_a… 0.134 816. -402. 12.2 8.29 -0.0171
## 2 Pocillopora_a… 0.134 816. -402. 12.2 8.29 -0.0171
## 3 Pocillopora_a… 0.134 816. -402. 12.2 8.29 -0.0171
## 4 Pocillopora_a… 0.134 816. -402. 12.2 8.29 -0.0171
## 5 Pocillopora_a… 0.134 816. -402. 12.2 8.29 -0.0171
## 6 Pocillopora_a… 0.134 816. -402. 12.2 8.29 -0.0171
## 7 Pocillopora_a… 0.134 816. -402. 12.2 8.29 -0.0171
## 8 Pocillopora_a… 0.134 816. -402. 12.2 8.29 -0.0171
## 9 Pocillopora_a… 0.134 816. -402. 12.2 8.29 -0.0171
## 10 Pocillopora_a… 0.134 816. -402. 12.2 8.29 -0.0171
## # ℹ 36 more rows
## # ℹ 29 more variables: OriginFlat <dbl>, TreatmentVariable.OriginFlat <dbl>,
## # se_.Intercept. <dbl>, se_TreatmentVariable <dbl>, se_OriginFlat <dbl>,
## # se_TreatmentVariable.OriginFlat <dbl>, Chisq_Treatment <dbl>,
## # Chisq_Origin <dbl>, Chisq_Treatment.Origin <dbl>, Df_Treatment <int>,
## # Df_Origin <int>, Df_Treatment.Origin <int>, P_Treatment <dbl>,
## # P_Origin <dbl>, P_Treatment.Origin <dbl>, Treatment <fct>, …library(nlme)
library(emmeans)
g25351.lme <- lme(Expression~Origin*Treatment2, random = ~1|Colony, data=biomin_all_filtered_long_g25351 , na.action=na.exclude)
car::Anova(g25351.lme, type=3)## Analysis of Deviance Table (Type III tests)
##
## Response: Expression
## Chisq Df Pr(>Chisq)
## (Intercept) 151.907 1 < 2.2e-16 ***
## Origin 10.058 1 0.001517 **
## Treatment2 1.959 1 0.161619
## Origin:Treatment2 1.039 1 0.308064
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1tukey3<- emmeans(g25351.lme, list(pairwise ~ Origin), adjust = "tukey")## NOTE: Results may be misleading due to involvement in interactionstukey3## $`emmeans of Origin`
## Origin emmean SE df lower.CL upper.CL
## RF 6454 354 19 5713 7195
## RS 3869 339 19 3159 4578
##
## Results are averaged over the levels of: Treatment2
## Degrees-of-freedom method: containment
## Confidence level used: 0.95
##
## $`pairwise differences of Origin`
## 1 estimate SE df t.ratio p.value
## RF - RS 2585 482 23 5.359 <.0001
##
## Results are averaged over the levels of: Treatment2
## Degrees-of-freedom method: containmentlibrary(Rmisc)
g25351_sum<-summarySE(biomin_all_filtered_long_g25351 , measurevar='Expression', groupvars=c('Origin', 'Treatment2'), na.rm=TRUE, conf.interval = 0.95)
g25351_sum## Origin Treatment2 N Expression sd se ci
## 1 RF Stable 11 5958.631 1717.4980 517.8451 1153.8309
## 2 RF Variable 11 6890.207 2342.7581 706.3681 1573.8863
## 3 RS Stable 12 3894.293 756.1492 218.2815 480.4342
## 4 RS Variable 12 3886.639 1300.8557 375.5247 826.5242Figure
pd<- position_dodge(0.2)
g25351_fig<-ggplot(data=g25351_sum, aes(y=Expression, x=Treatment2, color=Origin),group = interaction(Origin))+
geom_point(data=biomin_all_filtered_long_g25351,aes(y=Expression, x=Treatment2, color=Origin), alpha=0.4, position = pd)+
geom_line(aes(group = interaction(Origin), stat="identity"),position=position_dodge(0.2))+
geom_point(size=3, stat="identity", position = pd)+
geom_errorbar(aes(ymin=Expression-se, ymax=Expression+se), stat="identity",width=0.2, position = pd)+
scale_fill_manual("Origin", values=c("RS"='#26519B', "RF"= "#FE180C"))+
scale_color_manual("Origin", values=c("RS"='#26519B', "RF"= "#FE180C"))+
scale_y_continuous(expression(Mammalian~ependymin~related~protein~1-like))+
ggtitle(~MEBrown)+
theme_classic()+
theme(axis.text.x=element_text(vjust=0.5,size=12),#angling the labels on the x-axis
plot.title = element_text(margin = margin(t = 10, b = 10), hjust=0.5),#telling it where to position our plot title
panel.background= element_rect(fill=NA, color='black'),#this is making the black box around the graph
#strip.background = element_blank(),
#strip.text = element_blank(),
legend.title = element_text(vjust=0.5,size=12),
legend.position="none",
axis.text.y = element_text(vjust=0.5, size=12), #making the axis text larger
axis.title.x = element_blank(),#making the axis title larger
axis.title.y = element_text(size=12))#making the axis title larger ## Warning in geom_line(aes(group = interaction(Origin), stat = "identity"), :
## Ignoring unknown aesthetics: statg25351_fig
vitellogenin-like protein
biomin_all_filtered_long_g13222 <- biomin_all_filtered_long %>%
filter(Gene == "Pocillopora_acuta_HIv2___TS.g13222.t1b")
biomin_all_filtered_long_g13222 ## # A tibble: 184 × 36
## Gene Dispersion AIC logLik meanExp X.Intercept. TreatmentVariable
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Pocillopora_a… 0.0506 621. -305. 9.50 6.53 -0.123
## 2 Pocillopora_a… 0.0506 621. -305. 9.50 6.53 -0.123
## 3 Pocillopora_a… 0.0506 621. -305. 9.50 6.53 -0.123
## 4 Pocillopora_a… 0.0506 621. -305. 9.50 6.53 -0.123
## 5 Pocillopora_a… 0.0506 621. -305. 9.50 6.53 -0.123
## 6 Pocillopora_a… 0.0506 621. -305. 9.50 6.53 -0.123
## 7 Pocillopora_a… 0.0506 621. -305. 9.50 6.53 -0.123
## 8 Pocillopora_a… 0.0506 621. -305. 9.50 6.53 -0.123
## 9 Pocillopora_a… 0.0506 621. -305. 9.50 6.53 -0.123
## 10 Pocillopora_a… 0.0506 621. -305. 9.50 6.53 -0.123
## # ℹ 174 more rows
## # ℹ 29 more variables: OriginFlat <dbl>, TreatmentVariable.OriginFlat <dbl>,
## # se_.Intercept. <dbl>, se_TreatmentVariable <dbl>, se_OriginFlat <dbl>,
## # se_TreatmentVariable.OriginFlat <dbl>, Chisq_Treatment <dbl>,
## # Chisq_Origin <dbl>, Chisq_Treatment.Origin <dbl>, Df_Treatment <int>,
## # Df_Origin <int>, Df_Treatment.Origin <int>, P_Treatment <dbl>,
## # P_Origin <dbl>, P_Treatment.Origin <dbl>, Treatment <fct>, …library(nlme)
library(emmeans)
g13222.lme <- lme(Expression~Origin*Treatment2, random = ~1|Colony, data=biomin_all_filtered_long_g13222 , na.action=na.exclude)
car::Anova(g13222.lme, type=3)## Analysis of Deviance Table (Type III tests)
##
## Response: Expression
## Chisq Df Pr(>Chisq)
## (Intercept) 1126.7286 1 < 2.2e-16 ***
## Origin 41.2035 1 1.372e-10 ***
## Treatment2 2.9349 1 0.08669 .
## Origin:Treatment2 0.5477 1 0.45926
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1tukey3<- emmeans(g13222.lme, list(pairwise ~ Origin), adjust = "tukey")## NOTE: Results may be misleading due to involvement in interactionstukey3## $`emmeans of Origin`
## Origin emmean SE df lower.CL upper.CL
## RF 832 22.8 19 784 879
## RS 636 22.2 19 590 683
##
## Results are averaged over the levels of: Treatment2
## Degrees-of-freedom method: containment
## Confidence level used: 0.95
##
## $`pairwise differences of Origin`
## 1 estimate SE df t.ratio p.value
## RF - RS 195 24.2 161 8.053 <.0001
##
## Results are averaged over the levels of: Treatment2
## Degrees-of-freedom method: containmentlibrary(Rmisc)
g13222_sum<-summarySE(biomin_all_filtered_long_g13222 , measurevar='Expression', groupvars=c('Origin', 'Treatment2'), na.rm=TRUE, conf.interval = 0.95)
g13222_sum## Origin Treatment2 N Expression sd se ci
## 1 RF Stable 44 859.5886 158.58611 23.90776 48.21459
## 2 RF Variable 44 821.8027 114.92648 17.32582 34.94084
## 3 RS Stable 48 660.1489 146.62900 21.16407 42.57662
## 4 RS Variable 48 599.7645 70.60441 10.19087 20.50138Figure
pd<- position_dodge(0.2)
g13222_fig<-ggplot(data=g13222_sum, aes(y=Expression, x=Treatment2, color=Origin),group = interaction(Origin))+
geom_point(data=biomin_all_filtered_long_g13222,aes(y=Expression, x=Treatment2, color=Origin), alpha=0.4, position = pd)+
geom_line(aes(group = interaction(Origin), stat="identity"),position=position_dodge(0.2))+
geom_point(size=3, stat="identity", position = pd)+
geom_errorbar(aes(ymin=Expression-se, ymax=Expression+se), stat="identity",width=0.2, position = pd)+
scale_fill_manual("Origin", values=c("RS"='#26519B', "RF"= "#FE180C"))+
scale_color_manual("Origin", values=c("RS"='#26519B', "RF"= "#FE180C"))+
scale_y_continuous(expression(vitellogenin-like~~protein))+
ggtitle(~MEBrown)+
theme_classic()+
theme(axis.text.x=element_text(vjust=0.5,size=12),#angling the labels on the x-axis
plot.title = element_text(margin = margin(t = 10, b = 10), hjust=0.5),#telling it where to position our plot title
panel.background= element_rect(fill=NA, color='black'),#this is making the black box around the graph
#strip.background = element_blank(),
#strip.text = element_blank(),
legend.title = element_text(vjust=0.5,size=12),
legend.position="none",
axis.text.y = element_text(vjust=0.5, size=12), #making the axis text larger
axis.title.x = element_blank(),#making the axis title larger
axis.title.y = element_text(size=12))#making the axis title larger ## Warning in geom_line(aes(group = interaction(Origin), stat = "identity"), :
## Ignoring unknown aesthetics: statg13222_fig
uncharacterized protein LOC111323869
biomin_all_filtered_long_g21232 <- biomin_all_filtered_long %>%
filter(Gene == "Pocillopora_acuta_HIv2___RNAseq.g21232.t1")
biomin_all_filtered_long_g21232 ## # A tibble: 46 × 36
## Gene Dispersion AIC logLik meanExp X.Intercept. TreatmentVariable
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Pocillopora_a… 0.0605 353. -171. 5.22 3.54 -0.0600
## 2 Pocillopora_a… 0.0605 353. -171. 5.22 3.54 -0.0600
## 3 Pocillopora_a… 0.0605 353. -171. 5.22 3.54 -0.0600
## 4 Pocillopora_a… 0.0605 353. -171. 5.22 3.54 -0.0600
## 5 Pocillopora_a… 0.0605 353. -171. 5.22 3.54 -0.0600
## 6 Pocillopora_a… 0.0605 353. -171. 5.22 3.54 -0.0600
## 7 Pocillopora_a… 0.0605 353. -171. 5.22 3.54 -0.0600
## 8 Pocillopora_a… 0.0605 353. -171. 5.22 3.54 -0.0600
## 9 Pocillopora_a… 0.0605 353. -171. 5.22 3.54 -0.0600
## 10 Pocillopora_a… 0.0605 353. -171. 5.22 3.54 -0.0600
## # ℹ 36 more rows
## # ℹ 29 more variables: OriginFlat <dbl>, TreatmentVariable.OriginFlat <dbl>,
## # se_.Intercept. <dbl>, se_TreatmentVariable <dbl>, se_OriginFlat <dbl>,
## # se_TreatmentVariable.OriginFlat <dbl>, Chisq_Treatment <dbl>,
## # Chisq_Origin <dbl>, Chisq_Treatment.Origin <dbl>, Df_Treatment <int>,
## # Df_Origin <int>, Df_Treatment.Origin <int>, P_Treatment <dbl>,
## # P_Origin <dbl>, P_Treatment.Origin <dbl>, Treatment <fct>, …library(nlme)
library(emmeans)
g21232.lme <- lme(Expression~Origin*Treatment2, random = ~1|Colony, data=biomin_all_filtered_long_g21232 , na.action=na.exclude)
car::Anova(g21232.lme, type=3)## Analysis of Deviance Table (Type III tests)
##
## Response: Expression
## Chisq Df Pr(>Chisq)
## (Intercept) 209.9697 1 < 2e-16 ***
## Origin 4.7421 1 0.02943 *
## Treatment2 0.1307 1 0.71769
## Origin:Treatment2 0.0003 1 0.98700
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1tukey3<- emmeans(g21232.lme, list(pairwise ~ Origin), adjust = "tukey")## NOTE: Results may be misleading due to involvement in interactionstukey3## $`emmeans of Origin`
## Origin emmean SE df lower.CL upper.CL
## RF 42.0 2.4 19 36.9 47.0
## RS 33.7 2.3 19 28.8 38.5
##
## Results are averaged over the levels of: Treatment2
## Degrees-of-freedom method: containment
## Confidence level used: 0.95
##
## $`pairwise differences of Origin`
## 1 estimate SE df t.ratio p.value
## RF - RS 8.29 3.02 23 2.747 0.0115
##
## Results are averaged over the levels of: Treatment2
## Degrees-of-freedom method: containmentlibrary(Rmisc)
g21232_sum<-summarySE(biomin_all_filtered_long_g21232 , measurevar='Expression', groupvars=c('Origin', 'Treatment2'), na.rm=TRUE, conf.interval = 0.95)
g21232_sum## Origin Treatment2 N Expression sd se ci
## 1 RF Stable 11 42.29276 10.185331 3.070993 6.842599
## 2 RF Variable 11 41.06536 11.600467 3.497672 7.793300
## 3 RS Stable 12 33.84473 7.075632 2.042559 4.495642
## 4 RS Variable 12 32.69394 10.921097 3.152649 6.938934Figure
pd<- position_dodge(0.2)
g21232_fig<-ggplot(data=g21232_sum, aes(y=Expression, x=Treatment2, color=Origin),group = interaction(Origin))+
geom_point(data=biomin_all_filtered_long_g21232,aes(y=Expression, x=Treatment2, color=Origin), alpha=0.4, position = pd)+
geom_line(aes(group = interaction(Origin), stat="identity"),position=position_dodge(0.2))+
geom_point(size=3, stat="identity", position = pd)+
geom_errorbar(aes(ymin=Expression-se, ymax=Expression+se), stat="identity",width=0.2, position = pd)+
scale_fill_manual("Origin", values=c("RS"='#26519B', "RF"= "#FE180C"))+
scale_color_manual("Origin", values=c("RS"='#26519B', "RF"= "#FE180C"))+
scale_y_continuous(expression(uncharacterized~protein~LOC111323869))+
ggtitle(~MEBrown)+
theme_classic()+
theme(axis.text.x=element_text(vjust=0.5,size=12),#angling the labels on the x-axis
plot.title = element_text(margin = margin(t = 10, b = 10), hjust=0.5),#telling it where to position our plot title
panel.background= element_rect(fill=NA, color='black'),#this is making the black box around the graph
#strip.background = element_blank(),
#strip.text = element_blank(),
legend.title = element_text(vjust=0.5,size=12),
legend.position="none",
axis.text.y = element_text(vjust=0.5, size=12), #making the axis text larger
axis.title.x = element_blank(),#making the axis title larger
axis.title.y = element_text(size=12))#making the axis title larger ## Warning in geom_line(aes(group = interaction(Origin), stat = "identity"), :
## Ignoring unknown aesthetics: statg21232_fig
Cephalotoxin-like protein
biomin_all_filtered_long_g5013 <- biomin_all_filtered_long %>%
filter(Gene == "Pocillopora_acuta_HIv2___RNAseq.g5013.t1")
biomin_all_filtered_long_g5013 ## # A tibble: 46 × 36
## Gene Dispersion AIC logLik meanExp X.Intercept. TreatmentVariable
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Pocillopora_a… 0.471 304. -146. 3.16 1.83 0.206
## 2 Pocillopora_a… 0.471 304. -146. 3.16 1.83 0.206
## 3 Pocillopora_a… 0.471 304. -146. 3.16 1.83 0.206
## 4 Pocillopora_a… 0.471 304. -146. 3.16 1.83 0.206
## 5 Pocillopora_a… 0.471 304. -146. 3.16 1.83 0.206
## 6 Pocillopora_a… 0.471 304. -146. 3.16 1.83 0.206
## 7 Pocillopora_a… 0.471 304. -146. 3.16 1.83 0.206
## 8 Pocillopora_a… 0.471 304. -146. 3.16 1.83 0.206
## 9 Pocillopora_a… 0.471 304. -146. 3.16 1.83 0.206
## 10 Pocillopora_a… 0.471 304. -146. 3.16 1.83 0.206
## # ℹ 36 more rows
## # ℹ 29 more variables: OriginFlat <dbl>, TreatmentVariable.OriginFlat <dbl>,
## # se_.Intercept. <dbl>, se_TreatmentVariable <dbl>, se_OriginFlat <dbl>,
## # se_TreatmentVariable.OriginFlat <dbl>, Chisq_Treatment <dbl>,
## # Chisq_Origin <dbl>, Chisq_Treatment.Origin <dbl>, Df_Treatment <int>,
## # Df_Origin <int>, Df_Treatment.Origin <int>, P_Treatment <dbl>,
## # P_Origin <dbl>, P_Treatment.Origin <dbl>, Treatment <fct>, …library(nlme)
library(emmeans)
g5013.lme <- lme(Expression~Origin*Treatment2, random = ~1|Colony, data=biomin_all_filtered_long_g5013 , na.action=na.exclude)
car::Anova(g5013.lme, type=3)## Analysis of Deviance Table (Type III tests)
##
## Response: Expression
## Chisq Df Pr(>Chisq)
## (Intercept) 45.2530 1 1.732e-11 ***
## Origin 6.6236 1 0.01006 *
## Treatment2 0.0840 1 0.77190
## Origin:Treatment2 0.0561 1 0.81279
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1tukey3<- emmeans(g5013.lme, list(pairwise ~ Origin), adjust = "tukey")## NOTE: Results may be misleading due to involvement in interactionstukey3## $`emmeans of Origin`
## Origin emmean SE df lower.CL upper.CL
## RF 13.55 1.49 19 10.43 16.7
## RS 7.14 1.43 19 4.16 10.1
##
## Results are averaged over the levels of: Treatment2
## Degrees-of-freedom method: containment
## Confidence level used: 0.95
##
## $`pairwise differences of Origin`
## 1 estimate SE df t.ratio p.value
## RF - RS 6.41 1.98 23 3.237 0.0036
##
## Results are averaged over the levels of: Treatment2
## Degrees-of-freedom method: containmentlibrary(Rmisc)
g5013_sum<-summarySE(biomin_all_filtered_long_g5013 , measurevar='Expression', groupvars=c('Origin', 'Treatment2'), na.rm=TRUE, conf.interval = 0.95)
g5013_sum## Origin Treatment2 N Expression sd se ci
## 1 RF Stable 11 13.339697 7.878658 2.375505 5.292955
## 2 RF Variable 11 14.078906 6.431123 1.939056 4.320487
## 3 RS Stable 12 6.188827 5.219635 1.506779 3.316398
## 4 RS Variable 12 7.764125 6.413793 1.851503 4.075130Figure
pd<- position_dodge(0.2)
g5013_fig<-ggplot(data=g5013_sum, aes(y=Expression, x=Treatment2, color=Origin),group = interaction(Origin))+
geom_point(data=biomin_all_filtered_long_g5013,aes(y=Expression, x=Treatment2, color=Origin), alpha=0.4, position = pd)+
geom_line(aes(group = interaction(Origin), stat="identity"),position=position_dodge(0.2))+
geom_point(size=3, stat="identity", position = pd)+
geom_errorbar(aes(ymin=Expression-se, ymax=Expression+se), stat="identity",width=0.2, position = pd)+
scale_fill_manual("Origin", values=c("RS"='#26519B', "RF"= "#FE180C"))+
scale_color_manual("Origin", values=c("RS"='#26519B', "RF"= "#FE180C"))+
scale_y_continuous(expression(Cephalotoxin-like~protein))+
ggtitle(~MEBrown)+
theme_classic()+
theme(axis.text.x=element_text(vjust=0.5,size=12),#angling the labels on the x-axis
plot.title = element_text(margin = margin(t = 10, b = 10), hjust=0.5),#telling it where to position our plot title
panel.background= element_rect(fill=NA, color='black'),#this is making the black box around the graph
#strip.background = element_blank(),
#strip.text = element_blank(),
legend.title = element_text(vjust=0.5,size=12),
legend.position="none",
axis.text.y = element_text(vjust=0.5, size=12), #making the axis text larger
axis.title.x = element_blank(),#making the axis title larger
axis.title.y = element_text(size=12))#making the axis title larger ## Warning in geom_line(aes(group = interaction(Origin), stat = "identity"), :
## Ignoring unknown aesthetics: statg5013_fig
thioredoxin reductase 1, cytoplasmic-like
biomin_all_filtered_long_g10093 <- biomin_all_filtered_long %>%
filter(Gene == "Pocillopora_acuta_HIv2___RNAseq.g10093.t2")
biomin_all_filtered_long_g10093 ## # A tibble: 46 × 36
## Gene Dispersion AIC logLik meanExp X.Intercept. TreatmentVariable
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Pocillopora_a… 0.0876 385. -187. 5.74 3.88 -0.142
## 2 Pocillopora_a… 0.0876 385. -187. 5.74 3.88 -0.142
## 3 Pocillopora_a… 0.0876 385. -187. 5.74 3.88 -0.142
## 4 Pocillopora_a… 0.0876 385. -187. 5.74 3.88 -0.142
## 5 Pocillopora_a… 0.0876 385. -187. 5.74 3.88 -0.142
## 6 Pocillopora_a… 0.0876 385. -187. 5.74 3.88 -0.142
## 7 Pocillopora_a… 0.0876 385. -187. 5.74 3.88 -0.142
## 8 Pocillopora_a… 0.0876 385. -187. 5.74 3.88 -0.142
## 9 Pocillopora_a… 0.0876 385. -187. 5.74 3.88 -0.142
## 10 Pocillopora_a… 0.0876 385. -187. 5.74 3.88 -0.142
## # ℹ 36 more rows
## # ℹ 29 more variables: OriginFlat <dbl>, TreatmentVariable.OriginFlat <dbl>,
## # se_.Intercept. <dbl>, se_TreatmentVariable <dbl>, se_OriginFlat <dbl>,
## # se_TreatmentVariable.OriginFlat <dbl>, Chisq_Treatment <dbl>,
## # Chisq_Origin <dbl>, Chisq_Treatment.Origin <dbl>, Df_Treatment <int>,
## # Df_Origin <int>, Df_Treatment.Origin <int>, P_Treatment <dbl>,
## # P_Origin <dbl>, P_Treatment.Origin <dbl>, Treatment <fct>, …library(nlme)
library(emmeans)
g10093.lme <- lme(Expression~Origin*Treatment2, random = ~1|Colony, data=biomin_all_filtered_long_g10093 , na.action=na.exclude)
car::Anova(g10093.lme, type=3)## Analysis of Deviance Table (Type III tests)
##
## Response: Expression
## Chisq Df Pr(>Chisq)
## (Intercept) 205.6032 1 < 2.2e-16 ***
## Origin 10.8768 1 0.0009738 ***
## Treatment2 0.2141 1 0.6435866
## Origin:Treatment2 1.3642 1 0.2428065
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1tukey3<- emmeans(g10093.lme, list(pairwise ~ Origin), adjust = "tukey")## NOTE: Results may be misleading due to involvement in interactionstukey3## $`emmeans of Origin`
## Origin emmean SE df lower.CL upper.CL
## RF 66.4 3.89 19 58.2 74.5
## RS 43.9 3.76 19 36.0 51.7
##
## Results are averaged over the levels of: Treatment2
## Degrees-of-freedom method: containment
## Confidence level used: 0.95
##
## $`pairwise differences of Origin`
## 1 estimate SE df t.ratio p.value
## RF - RS 22.5 4.63 23 4.862 0.0001
##
## Results are averaged over the levels of: Treatment2
## Degrees-of-freedom method: containmentlibrary(Rmisc)
g10093_sum<-summarySE(biomin_all_filtered_long_g10093 , measurevar='Expression', groupvars=c('Origin', 'Treatment2'), na.rm=TRUE, conf.interval = 0.95)
g10093_sum## Origin Treatment2 N Expression sd se ci
## 1 RF Stable 11 63.79961 16.37758 4.938025 11.002606
## 2 RF Variable 11 65.98269 22.07239 6.655076 14.828433
## 3 RS Stable 12 47.40346 10.98308 3.170542 6.978315
## 4 RS Variable 12 41.95704 10.53729 3.041854 6.695076Figure
pd<- position_dodge(0.2)
g10093_fig<-ggplot(data=g10093_sum, aes(y=Expression, x=Treatment2, color=Origin),group = interaction(Origin))+
geom_point(data=biomin_all_filtered_long_g10093,aes(y=Expression, x=Treatment2, color=Origin), alpha=0.4, position = pd)+
geom_line(aes(group = interaction(Origin), stat="identity"),position=position_dodge(0.2))+
geom_point(size=3, stat="identity", position = pd)+
geom_errorbar(aes(ymin=Expression-se, ymax=Expression+se), stat="identity",width=0.2, position = pd)+
scale_fill_manual("Origin", values=c("RS"='#26519B', "RF"= "#FE180C"))+
scale_color_manual("Origin", values=c("RS"='#26519B', "RF"= "#FE180C"))+
scale_y_continuous(expression(thioredoxin~reductase~1))+
ggtitle(~MEBrown)+
theme_classic()+
theme(axis.text.x=element_text(vjust=0.5,size=12),#angling the labels on the x-axis
plot.title = element_text(margin = margin(t = 10, b = 10), hjust=0.5),#telling it where to position our plot title
panel.background= element_rect(fill=NA, color='black'),#this is making the black box around the graph
#strip.background = element_blank(),
#strip.text = element_blank(),
legend.title = element_text(vjust=0.5,size=12),
legend.position="none",
axis.text.y = element_text(vjust=0.5, size=12), #making the axis text larger
axis.title.x = element_blank(),#making the axis title larger
axis.title.y = element_text(size=12))#making the axis title larger ## Warning in geom_line(aes(group = interaction(Origin), stat = "identity"), :
## Ignoring unknown aesthetics: statg10093_fig
carbonic anhydrase (STPCA2-2)
biomin_all_filtered_long_g13824 <- biomin_all_filtered_long %>%
filter(Gene == "Pocillopora_acuta_HIv2___RNAseq.g13824.t1")
biomin_all_filtered_long_g13824 ## # A tibble: 46 × 36
## Gene Dispersion AIC logLik meanExp X.Intercept. TreatmentVariable
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Pocillopora_a… 2.35 453. -220. 4.66 2.36 -0.0648
## 2 Pocillopora_a… 2.35 453. -220. 4.66 2.36 -0.0648
## 3 Pocillopora_a… 2.35 453. -220. 4.66 2.36 -0.0648
## 4 Pocillopora_a… 2.35 453. -220. 4.66 2.36 -0.0648
## 5 Pocillopora_a… 2.35 453. -220. 4.66 2.36 -0.0648
## 6 Pocillopora_a… 2.35 453. -220. 4.66 2.36 -0.0648
## 7 Pocillopora_a… 2.35 453. -220. 4.66 2.36 -0.0648
## 8 Pocillopora_a… 2.35 453. -220. 4.66 2.36 -0.0648
## 9 Pocillopora_a… 2.35 453. -220. 4.66 2.36 -0.0648
## 10 Pocillopora_a… 2.35 453. -220. 4.66 2.36 -0.0648
## # ℹ 36 more rows
## # ℹ 29 more variables: OriginFlat <dbl>, TreatmentVariable.OriginFlat <dbl>,
## # se_.Intercept. <dbl>, se_TreatmentVariable <dbl>, se_OriginFlat <dbl>,
## # se_TreatmentVariable.OriginFlat <dbl>, Chisq_Treatment <dbl>,
## # Chisq_Origin <dbl>, Chisq_Treatment.Origin <dbl>, Df_Treatment <int>,
## # Df_Origin <int>, Df_Treatment.Origin <int>, P_Treatment <dbl>,
## # P_Origin <dbl>, P_Treatment.Origin <dbl>, Treatment <fct>, …library(nlme)
library(emmeans)
g13824.lme <- lme(Expression~Origin*Treatment2, random = ~1|Colony, data=biomin_all_filtered_long_g13824 , na.action=na.exclude)
car::Anova(g13824.lme, type=3)## Analysis of Deviance Table (Type III tests)
##
## Response: Expression
## Chisq Df Pr(>Chisq)
## (Intercept) 77.9228 1 <2e-16 ***
## Origin 78.0753 1 <2e-16 ***
## Treatment2 2.0231 1 0.1549
## Origin:Treatment2 1.0392 1 0.3080
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1tukey3<- emmeans(g13824.lme, list(pairwise ~ Origin), adjust = "tukey")## NOTE: Results may be misleading due to involvement in interactionstukey3## $`emmeans of Origin`
## Origin emmean SE df lower.CL upper.CL
## RF 158.76 17.7 19 121.7 195.8
## RS -5.68 17.3 19 -41.9 30.5
##
## Results are averaged over the levels of: Treatment2
## Degrees-of-freedom method: containment
## Confidence level used: 0.95
##
## $`pairwise differences of Origin`
## 1 estimate SE df t.ratio p.value
## RF - RS 164 17 23 9.671 <.0001
##
## Results are averaged over the levels of: Treatment2
## Degrees-of-freedom method: containmentlibrary(Rmisc)
g13824_sum<-summarySE(biomin_all_filtered_long_g13824 , measurevar='Expression', groupvars=c('Origin', 'Treatment2'), na.rm=TRUE, conf.interval = 0.95)
g13824_sum## Origin Treatment2 N Expression sd se ci
## 1 RF Stable 11 156.36191 109.769844 33.096853 73.744384
## 2 RF Variable 11 135.61447 87.446714 26.366176 58.747502
## 3 RS Stable 12 10.39836 9.661042 2.788902 6.138333
## 4 RS Variable 12 10.23764 8.081504 2.332929 5.134743Figure
pd<- position_dodge(0.2)
g13824_fig<-ggplot(data=g13824_sum, aes(y=Expression, x=Treatment2, color=Origin),group = interaction(Origin))+
geom_point(data=biomin_all_filtered_long_g13824,aes(y=Expression, x=Treatment2, color=Origin), alpha=0.4, position = pd)+
geom_line(aes(group = interaction(Origin), stat="identity"),position=position_dodge(0.2))+
geom_point(size=3, stat="identity", position = pd)+
geom_errorbar(aes(ymin=Expression-se, ymax=Expression+se), stat="identity",width=0.2, position = pd)+
scale_fill_manual("Origin", values=c("RS"='#26519B', "RF"= "#FE180C"))+
scale_color_manual("Origin", values=c("RS"='#26519B', "RF"= "#FE180C"))+
scale_y_continuous(expression(carbonic~anhydrase~(STPCA2-2)))+
ggtitle(~MEBrown)+
theme_classic()+
theme(axis.text.x=element_text(vjust=0.5,size=12),#angling the labels on the x-axis
plot.title = element_text(margin = margin(t = 10, b = 10), hjust=0.5),#telling it where to position our plot title
panel.background= element_rect(fill=NA, color='black'),#this is making the black box around the graph
#strip.background = element_blank(),
#strip.text = element_blank(),
legend.title = element_text(vjust=0.5,size=12),
legend.position="none",
axis.text.y = element_text(vjust=0.5, size=12), #making the axis text larger
axis.title.x = element_blank(),#making the axis title larger
axis.title.y = element_text(size=12))#making the axis title larger ## Warning in geom_line(aes(group = interaction(Origin), stat = "identity"), :
## Ignoring unknown aesthetics: statg13824_fig
Comparison of all Brown and Red sig. genes
biomin_compare_figs <- cowplot::plot_grid(g12304_fig,g7908_fig,g25351_fig,g13222_fig,g21232_fig, g5013_fig,g10093_fig,g13824_fig,nrow=3)
biomin_compare_figs
ggsave(filename="../../output/WGCNA/GO_analysis/biomin_Brown_Red_compare_figs.png", plot=biomin_compare_figs, dpi=300, height=8, units="in", limitsize=FALSE)## Saving 7 x 8 in imageOld Pink and Blue code below:
MEBlue signficant genes
Thioredoxin reductase 1
biomin_all_filtered_long_g10093 <- biomin_all_filtered_long %>%
filter(Gene == "Pocillopora_acuta_HIv2___RNAseq.g10093.t2")
biomin_all_filtered_long_g10093 ## # A tibble: 46 × 36
## Gene Dispersion AIC logLik meanExp X.Intercept. TreatmentVariable
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Pocillopora_a… 0.0876 385. -187. 5.74 3.88 -0.142
## 2 Pocillopora_a… 0.0876 385. -187. 5.74 3.88 -0.142
## 3 Pocillopora_a… 0.0876 385. -187. 5.74 3.88 -0.142
## 4 Pocillopora_a… 0.0876 385. -187. 5.74 3.88 -0.142
## 5 Pocillopora_a… 0.0876 385. -187. 5.74 3.88 -0.142
## 6 Pocillopora_a… 0.0876 385. -187. 5.74 3.88 -0.142
## 7 Pocillopora_a… 0.0876 385. -187. 5.74 3.88 -0.142
## 8 Pocillopora_a… 0.0876 385. -187. 5.74 3.88 -0.142
## 9 Pocillopora_a… 0.0876 385. -187. 5.74 3.88 -0.142
## 10 Pocillopora_a… 0.0876 385. -187. 5.74 3.88 -0.142
## # ℹ 36 more rows
## # ℹ 29 more variables: OriginFlat <dbl>, TreatmentVariable.OriginFlat <dbl>,
## # se_.Intercept. <dbl>, se_TreatmentVariable <dbl>, se_OriginFlat <dbl>,
## # se_TreatmentVariable.OriginFlat <dbl>, Chisq_Treatment <dbl>,
## # Chisq_Origin <dbl>, Chisq_Treatment.Origin <dbl>, Df_Treatment <int>,
## # Df_Origin <int>, Df_Treatment.Origin <int>, P_Treatment <dbl>,
## # P_Origin <dbl>, P_Treatment.Origin <dbl>, Treatment <fct>, …library(nlme)
library(emmeans)
g10093.lme <- lme(Expression~Origin*Treatment2, random = ~1|Colony, data=biomin_all_filtered_long_g10093 , na.action=na.exclude)
car::Anova(g10093.lme, type=3)## Analysis of Deviance Table (Type III tests)
##
## Response: Expression
## Chisq Df Pr(>Chisq)
## (Intercept) 205.6032 1 < 2.2e-16 ***
## Origin 10.8768 1 0.0009738 ***
## Treatment2 0.2141 1 0.6435866
## Origin:Treatment2 1.3642 1 0.2428065
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1tukey3<- emmeans(g10093.lme, list(pairwise ~ Origin), adjust = "tukey")## NOTE: Results may be misleading due to involvement in interactionstukey3## $`emmeans of Origin`
## Origin emmean SE df lower.CL upper.CL
## RF 66.4 3.89 19 58.2 74.5
## RS 43.9 3.76 19 36.0 51.7
##
## Results are averaged over the levels of: Treatment2
## Degrees-of-freedom method: containment
## Confidence level used: 0.95
##
## $`pairwise differences of Origin`
## 1 estimate SE df t.ratio p.value
## RF - RS 22.5 4.63 23 4.862 0.0001
##
## Results are averaged over the levels of: Treatment2
## Degrees-of-freedom method: containmentlibrary(Rmisc)
g10093_sum<-summarySE(biomin_all_filtered_long_g10093 , measurevar='Expression', groupvars=c('Origin', 'Treatment2'), na.rm=TRUE, conf.interval = 0.95)
g10093_sum## Origin Treatment2 N Expression sd se ci
## 1 RF Stable 11 63.79961 16.37758 4.938025 11.002606
## 2 RF Variable 11 65.98269 22.07239 6.655076 14.828433
## 3 RS Stable 12 47.40346 10.98308 3.170542 6.978315
## 4 RS Variable 12 41.95704 10.53729 3.041854 6.695076Figure
pd<- position_dodge(0.2)
g10093_fig<-ggplot(data=g10093_sum, aes(y=Expression, x=Treatment2, color=Origin),group = interaction(Origin))+
geom_point(data=biomin_all_filtered_long_g10093,aes(y=Expression, x=Treatment2, color=Origin), alpha=0.4, position = pd)+
geom_line(aes(group = interaction(Origin), stat="identity"),position=position_dodge(0.2))+
geom_point(size=3, stat="identity", position = pd)+
geom_errorbar(aes(ymin=Expression-se, ymax=Expression+se), stat="identity",width=0.2, position = pd)+
scale_fill_manual("Origin", values=c("RS"='#26519B', "RF"= "#FE180C"))+
scale_color_manual("Origin", values=c("RS"='#26519B', "RF"= "#FE180C"))+
scale_y_continuous(expression(Thioredoxin~reductase~1))+
ggtitle(~blue)+
theme_classic()+
theme(axis.text.x=element_text(vjust=0.5,size=12),#angling the labels on the x-axis
plot.title = element_text(margin = margin(t = 10, b = 10), hjust=0.5),#telling it where to position our plot title
panel.background= element_rect(fill=NA, color='black'),#this is making the black box around the graph
#strip.background = element_blank(),
#strip.text = element_blank(),
legend.title = element_text(vjust=0.5,size=12),
legend.position="none",
axis.text.y = element_text(vjust=0.5, size=12), #making the axis text larger
axis.title.x = element_blank(),#making the axis title larger
axis.title.y = element_text(size=12))#making the axis title larger ## Warning in geom_line(aes(group = interaction(Origin), stat = "identity"), :
## Ignoring unknown aesthetics: statg10093_fig
Carbonic anhydrase - STPCA2-2
biomin_all_filtered_long_g13824 <- biomin_all_filtered_long %>%
filter(Gene == "Pocillopora_acuta_HIv2___RNAseq.g13824.t1")
biomin_all_filtered_long_g13824 ## # A tibble: 46 × 36
## Gene Dispersion AIC logLik meanExp X.Intercept. TreatmentVariable
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Pocillopora_a… 2.35 453. -220. 4.66 2.36 -0.0648
## 2 Pocillopora_a… 2.35 453. -220. 4.66 2.36 -0.0648
## 3 Pocillopora_a… 2.35 453. -220. 4.66 2.36 -0.0648
## 4 Pocillopora_a… 2.35 453. -220. 4.66 2.36 -0.0648
## 5 Pocillopora_a… 2.35 453. -220. 4.66 2.36 -0.0648
## 6 Pocillopora_a… 2.35 453. -220. 4.66 2.36 -0.0648
## 7 Pocillopora_a… 2.35 453. -220. 4.66 2.36 -0.0648
## 8 Pocillopora_a… 2.35 453. -220. 4.66 2.36 -0.0648
## 9 Pocillopora_a… 2.35 453. -220. 4.66 2.36 -0.0648
## 10 Pocillopora_a… 2.35 453. -220. 4.66 2.36 -0.0648
## # ℹ 36 more rows
## # ℹ 29 more variables: OriginFlat <dbl>, TreatmentVariable.OriginFlat <dbl>,
## # se_.Intercept. <dbl>, se_TreatmentVariable <dbl>, se_OriginFlat <dbl>,
## # se_TreatmentVariable.OriginFlat <dbl>, Chisq_Treatment <dbl>,
## # Chisq_Origin <dbl>, Chisq_Treatment.Origin <dbl>, Df_Treatment <int>,
## # Df_Origin <int>, Df_Treatment.Origin <int>, P_Treatment <dbl>,
## # P_Origin <dbl>, P_Treatment.Origin <dbl>, Treatment <fct>, …library(nlme)
g13824.lme <- lme(Expression~Origin*Treatment2, random = ~1|Colony, data=biomin_all_filtered_long_g13824 , na.action=na.exclude)
car::Anova(g13824.lme, type=3)## Analysis of Deviance Table (Type III tests)
##
## Response: Expression
## Chisq Df Pr(>Chisq)
## (Intercept) 77.9228 1 <2e-16 ***
## Origin 78.0753 1 <2e-16 ***
## Treatment2 2.0231 1 0.1549
## Origin:Treatment2 1.0392 1 0.3080
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1tukey3<- emmeans(g13824.lme, list(pairwise ~ Origin), adjust = "tukey")## NOTE: Results may be misleading due to involvement in interactionstukey3## $`emmeans of Origin`
## Origin emmean SE df lower.CL upper.CL
## RF 158.76 17.7 19 121.7 195.8
## RS -5.68 17.3 19 -41.9 30.5
##
## Results are averaged over the levels of: Treatment2
## Degrees-of-freedom method: containment
## Confidence level used: 0.95
##
## $`pairwise differences of Origin`
## 1 estimate SE df t.ratio p.value
## RF - RS 164 17 23 9.671 <.0001
##
## Results are averaged over the levels of: Treatment2
## Degrees-of-freedom method: containmentlibrary(Rmisc)
g13824_sum<-summarySE(biomin_all_filtered_long_g13824 , measurevar='Expression', groupvars=c('Origin', 'Treatment2'), na.rm=TRUE, conf.interval = 0.95)
g13824_sum## Origin Treatment2 N Expression sd se ci
## 1 RF Stable 11 156.36191 109.769844 33.096853 73.744384
## 2 RF Variable 11 135.61447 87.446714 26.366176 58.747502
## 3 RS Stable 12 10.39836 9.661042 2.788902 6.138333
## 4 RS Variable 12 10.23764 8.081504 2.332929 5.134743###Figure
pd<- position_dodge(0.2)
g13824_fig<-ggplot(data=g13824_sum, aes(y=Expression, x=Treatment2, color=Origin),group = interaction(Origin))+
geom_point(data=biomin_all_filtered_long_g13824,aes(y=Expression, x=Treatment2, color=Origin), alpha=0.4, position = pd)+
geom_line(aes(group = interaction(Origin), stat="identity"),position=position_dodge(0.2))+
geom_point(size=3, stat="identity", position = pd)+
geom_errorbar(aes(ymin=Expression-se, ymax=Expression+se), stat="identity",width=0.2, position = pd)+
scale_fill_manual("Origin", values=c("RS"='#26519B', "RF"= "#FE180C"))+
scale_color_manual("Origin", values=c("RS"='#26519B', "RF"= "#FE180C"))+
scale_y_continuous(expression(Carbonic~anhydrase~("STPCA2-2")), limits=c(0,300))+
ggtitle(~blue)+
theme_classic()+
theme(axis.text.x=element_text(vjust=0.5,size=12),#angling the labels on the x-axis
plot.title = element_text(margin = margin(t = 10, b = 10), hjust=0.5),#telling it where to position our plot title
panel.background= element_rect(fill=NA, color='black'),#this is making the black box around the graph
#strip.background = element_blank(),
#strip.text = element_blank(),
legend.title = element_text(vjust=0.5,size=12),
legend.position="none",
axis.text.y = element_text(vjust=0.5, size=12), #making the axis text larger
axis.title.x = element_blank(),#making the axis title larger
axis.title.y = element_text(size=12))#making the axis title larger ## Warning in geom_line(aes(group = interaction(Origin), stat = "identity"), :
## Ignoring unknown aesthetics: statg13824_fig## Warning: Removed 1 rows containing missing values (`geom_point()`).
Mammalian ependymin-related protein 1-like
biomin_all_filtered_long_g25351 <- biomin_all_filtered_long %>%
filter(Gene == "Pocillopora_acuta_HIv2___RNAseq.g25351.t1")
biomin_all_filtered_long_g25351 ## # A tibble: 46 × 36
## Gene Dispersion AIC logLik meanExp X.Intercept. TreatmentVariable
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Pocillopora_a… 0.134 816. -402. 12.2 8.29 -0.0171
## 2 Pocillopora_a… 0.134 816. -402. 12.2 8.29 -0.0171
## 3 Pocillopora_a… 0.134 816. -402. 12.2 8.29 -0.0171
## 4 Pocillopora_a… 0.134 816. -402. 12.2 8.29 -0.0171
## 5 Pocillopora_a… 0.134 816. -402. 12.2 8.29 -0.0171
## 6 Pocillopora_a… 0.134 816. -402. 12.2 8.29 -0.0171
## 7 Pocillopora_a… 0.134 816. -402. 12.2 8.29 -0.0171
## 8 Pocillopora_a… 0.134 816. -402. 12.2 8.29 -0.0171
## 9 Pocillopora_a… 0.134 816. -402. 12.2 8.29 -0.0171
## 10 Pocillopora_a… 0.134 816. -402. 12.2 8.29 -0.0171
## # ℹ 36 more rows
## # ℹ 29 more variables: OriginFlat <dbl>, TreatmentVariable.OriginFlat <dbl>,
## # se_.Intercept. <dbl>, se_TreatmentVariable <dbl>, se_OriginFlat <dbl>,
## # se_TreatmentVariable.OriginFlat <dbl>, Chisq_Treatment <dbl>,
## # Chisq_Origin <dbl>, Chisq_Treatment.Origin <dbl>, Df_Treatment <int>,
## # Df_Origin <int>, Df_Treatment.Origin <int>, P_Treatment <dbl>,
## # P_Origin <dbl>, P_Treatment.Origin <dbl>, Treatment <fct>, …library(nlme)
g25351.lme <- lme(Expression~Origin*Treatment2, random = ~1|Colony, data=biomin_all_filtered_long_g25351 , na.action=na.exclude)
car::Anova(g25351.lme, type=3)## Analysis of Deviance Table (Type III tests)
##
## Response: Expression
## Chisq Df Pr(>Chisq)
## (Intercept) 151.907 1 < 2.2e-16 ***
## Origin 10.058 1 0.001517 **
## Treatment2 1.959 1 0.161619
## Origin:Treatment2 1.039 1 0.308064
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1tukey3<- emmeans(g25351.lme, list(pairwise ~ Origin), adjust = "tukey")## NOTE: Results may be misleading due to involvement in interactionstukey3## $`emmeans of Origin`
## Origin emmean SE df lower.CL upper.CL
## RF 6454 354 19 5713 7195
## RS 3869 339 19 3159 4578
##
## Results are averaged over the levels of: Treatment2
## Degrees-of-freedom method: containment
## Confidence level used: 0.95
##
## $`pairwise differences of Origin`
## 1 estimate SE df t.ratio p.value
## RF - RS 2585 482 23 5.359 <.0001
##
## Results are averaged over the levels of: Treatment2
## Degrees-of-freedom method: containmentlibrary(Rmisc)
g25351_sum<-summarySE(biomin_all_filtered_long_g25351 , measurevar='Expression', groupvars=c('Origin', 'Treatment2'), na.rm=TRUE, conf.interval = 0.95)
g25351_sum## Origin Treatment2 N Expression sd se ci
## 1 RF Stable 11 5958.631 1717.4980 517.8451 1153.8309
## 2 RF Variable 11 6890.207 2342.7581 706.3681 1573.8863
## 3 RS Stable 12 3894.293 756.1492 218.2815 480.4342
## 4 RS Variable 12 3886.639 1300.8557 375.5247 826.5242Figure
pd<- position_dodge(0.2)
g25351_fig<-ggplot(data=g25351_sum, aes(y=Expression, x=Treatment2, color=Origin),group = interaction(Origin))+
geom_point(data=biomin_all_filtered_long_g25351,aes(y=Expression, x=Treatment2, color=Origin), alpha=0.4, position = pd)+
geom_line(aes(group = interaction(Origin), stat="identity"),position=position_dodge(0.2))+
geom_point(size=3, stat="identity", position = pd)+
geom_errorbar(aes(ymin=Expression-se, ymax=Expression+se), stat="identity",width=0.2, position = pd)+
scale_fill_manual("Origin", values=c("RS"='#26519B', "RF"= "#FE180C"))+
scale_color_manual("Origin", values=c("RS"='#26519B', "RF"= "#FE180C"))+
scale_y_continuous(expression(Mammalian~ependymin-related~protein~1-like),limits=c(2000,10000))+
ggtitle(~blue)+
theme_classic()+
theme(axis.text.x=element_text(vjust=0.5,size=12),#angling the labels on the x-axis
plot.title = element_text(margin = margin(t = 10, b = 10), hjust=0.5),#telling it where to position our plot title
panel.background= element_rect(fill=NA, color='black'),#this is making the black box around the graph
#strip.background = element_blank(),
#strip.text = element_blank(),
legend.title = element_text(vjust=0.5,size=12),
legend.position="none",
axis.text.y = element_text(vjust=0.5, size=12), #making the axis text larger
axis.title.x = element_blank(),#making the axis title larger
axis.title.y = element_text(size=12))#making the axis title larger ## Warning in geom_line(aes(group = interaction(Origin), stat = "identity"), :
## Ignoring unknown aesthetics: statg25351_fig## Warning: Removed 1 rows containing missing values (`geom_point()`).
Cephalotoxin-like protein
biomin_all_filtered_long_g5013 <- biomin_all_filtered_long %>%
filter(Gene == "Pocillopora_acuta_HIv2___RNAseq.g5013.t1")
biomin_all_filtered_long_g5013 ## # A tibble: 46 × 36
## Gene Dispersion AIC logLik meanExp X.Intercept. TreatmentVariable
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Pocillopora_a… 0.471 304. -146. 3.16 1.83 0.206
## 2 Pocillopora_a… 0.471 304. -146. 3.16 1.83 0.206
## 3 Pocillopora_a… 0.471 304. -146. 3.16 1.83 0.206
## 4 Pocillopora_a… 0.471 304. -146. 3.16 1.83 0.206
## 5 Pocillopora_a… 0.471 304. -146. 3.16 1.83 0.206
## 6 Pocillopora_a… 0.471 304. -146. 3.16 1.83 0.206
## 7 Pocillopora_a… 0.471 304. -146. 3.16 1.83 0.206
## 8 Pocillopora_a… 0.471 304. -146. 3.16 1.83 0.206
## 9 Pocillopora_a… 0.471 304. -146. 3.16 1.83 0.206
## 10 Pocillopora_a… 0.471 304. -146. 3.16 1.83 0.206
## # ℹ 36 more rows
## # ℹ 29 more variables: OriginFlat <dbl>, TreatmentVariable.OriginFlat <dbl>,
## # se_.Intercept. <dbl>, se_TreatmentVariable <dbl>, se_OriginFlat <dbl>,
## # se_TreatmentVariable.OriginFlat <dbl>, Chisq_Treatment <dbl>,
## # Chisq_Origin <dbl>, Chisq_Treatment.Origin <dbl>, Df_Treatment <int>,
## # Df_Origin <int>, Df_Treatment.Origin <int>, P_Treatment <dbl>,
## # P_Origin <dbl>, P_Treatment.Origin <dbl>, Treatment <fct>, …library(nlme)
g5013.lme <- lme(Expression~Origin*Treatment2, random = ~1|Colony, data=biomin_all_filtered_long_g5013 , na.action=na.exclude)
car::Anova(g5013.lme, type=3)## Analysis of Deviance Table (Type III tests)
##
## Response: Expression
## Chisq Df Pr(>Chisq)
## (Intercept) 45.2530 1 1.732e-11 ***
## Origin 6.6236 1 0.01006 *
## Treatment2 0.0840 1 0.77190
## Origin:Treatment2 0.0561 1 0.81279
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1tukey3<- emmeans(g5013.lme, list(pairwise ~ Origin), adjust = "tukey")## NOTE: Results may be misleading due to involvement in interactionstukey3## $`emmeans of Origin`
## Origin emmean SE df lower.CL upper.CL
## RF 13.55 1.49 19 10.43 16.7
## RS 7.14 1.43 19 4.16 10.1
##
## Results are averaged over the levels of: Treatment2
## Degrees-of-freedom method: containment
## Confidence level used: 0.95
##
## $`pairwise differences of Origin`
## 1 estimate SE df t.ratio p.value
## RF - RS 6.41 1.98 23 3.237 0.0036
##
## Results are averaged over the levels of: Treatment2
## Degrees-of-freedom method: containmentlibrary(Rmisc)
g5013_sum<-summarySE(biomin_all_filtered_long_g5013 , measurevar='Expression', groupvars=c('Origin', 'Treatment2'), na.rm=TRUE, conf.interval = 0.95)
g5013_sum## Origin Treatment2 N Expression sd se ci
## 1 RF Stable 11 13.339697 7.878658 2.375505 5.292955
## 2 RF Variable 11 14.078906 6.431123 1.939056 4.320487
## 3 RS Stable 12 6.188827 5.219635 1.506779 3.316398
## 4 RS Variable 12 7.764125 6.413793 1.851503 4.075130Figure
pd<- position_dodge(0.2)
g5013_fig<-ggplot(data=g5013_sum, aes(y=Expression, x=Treatment2, color=Origin),group = interaction(Origin))+
geom_point(data=biomin_all_filtered_long_g5013,aes(y=Expression, x=Treatment2, color=Origin), alpha=0.4, position = pd)+
geom_line(aes(group = interaction(Origin), stat="identity"),position=position_dodge(0.2))+
geom_point(size=3, stat="identity", position = pd)+
geom_errorbar(aes(ymin=Expression-se, ymax=Expression+se), stat="identity",width=0.2, position = pd)+
scale_fill_manual("Origin", values=c("RS"='#26519B', "RF"= "#FE180C"))+
scale_color_manual("Origin", values=c("RS"='#26519B', "RF"= "#FE180C"))+
scale_y_continuous(expression(Cephalotoxin-like~protein), breaks=c(0,5,10,15,20,25))+
ggtitle(~blue)+
theme_classic()+
theme(axis.text.x=element_text(vjust=0.5,size=12),#angling the labels on the x-axis
plot.title = element_text(margin = margin(t = 10, b = 10), hjust=0.5),#telling it where to position our plot title
panel.background= element_rect(fill=NA, color='black'),#this is making the black box around the graph
#strip.background = element_blank(),
#strip.text = element_blank(),
legend.title = element_text(vjust=0.5,size=12),
legend.position="none",
axis.text.y = element_text(vjust=0.5, size=12), #making the axis text larger
axis.title.x = element_blank(),#making the axis title larger
axis.title.y = element_text(size=12))#making the axis title larger ## Warning in geom_line(aes(group = interaction(Origin), stat = "identity"), :
## Ignoring unknown aesthetics: statg5013_fig
MEPink signficant genes
Carbonic anhydrase - STPCA2-1
biomin_all_filtered_long_g12304 <- biomin_all_filtered_long %>%
filter(Gene == "Pocillopora_acuta_HIv2___TS.g12304.t1")
biomin_all_filtered_long_g12304 ## # A tibble: 184 × 36
## Gene Dispersion AIC logLik meanExp X.Intercept. TreatmentVariable
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Pocillopora_a… 0.457 841. -415. 11.6 7.88 0.0448
## 2 Pocillopora_a… 0.457 841. -415. 11.6 7.88 0.0448
## 3 Pocillopora_a… 0.457 841. -415. 11.6 7.88 0.0448
## 4 Pocillopora_a… 0.457 841. -415. 11.6 7.88 0.0448
## 5 Pocillopora_a… 0.457 841. -415. 11.6 7.88 0.0448
## 6 Pocillopora_a… 0.457 841. -415. 11.6 7.88 0.0448
## 7 Pocillopora_a… 0.457 841. -415. 11.6 7.88 0.0448
## 8 Pocillopora_a… 0.457 841. -415. 11.6 7.88 0.0448
## 9 Pocillopora_a… 0.457 841. -415. 11.6 7.88 0.0448
## 10 Pocillopora_a… 0.457 841. -415. 11.6 7.88 0.0448
## # ℹ 174 more rows
## # ℹ 29 more variables: OriginFlat <dbl>, TreatmentVariable.OriginFlat <dbl>,
## # se_.Intercept. <dbl>, se_TreatmentVariable <dbl>, se_OriginFlat <dbl>,
## # se_TreatmentVariable.OriginFlat <dbl>, Chisq_Treatment <dbl>,
## # Chisq_Origin <dbl>, Chisq_Treatment.Origin <dbl>, Df_Treatment <int>,
## # Df_Origin <int>, Df_Treatment.Origin <int>, P_Treatment <dbl>,
## # P_Origin <dbl>, P_Treatment.Origin <dbl>, Treatment <fct>, …library(nlme)
g12304.lme <- lme(Expression~Origin*Treatment2, random = ~1|Colony, data=biomin_all_filtered_long_g12304 , na.action=na.exclude)
car::Anova(g12304.lme, type=3)## Analysis of Deviance Table (Type III tests)
##
## Response: Expression
## Chisq Df Pr(>Chisq)
## (Intercept) 196.964 1 < 2.2e-16 ***
## Origin 55.461 1 9.531e-14 ***
## Treatment2 14.008 1 0.0001820 ***
## Origin:Treatment2 11.399 1 0.0007348 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1tukey3<- emmeans(g12304.lme, list(pairwise ~ Origin), adjust = "tukey")## NOTE: Results may be misleading due to involvement in interactionstukey3## $`emmeans of Origin`
## Origin emmean SE df lower.CL upper.CL
## RF 5238 384 19 4434 6042
## RS 2770 375 19 1985 3555
##
## Results are averaged over the levels of: Treatment2
## Degrees-of-freedom method: containment
## Confidence level used: 0.95
##
## $`pairwise differences of Origin`
## 1 estimate SE df t.ratio p.value
## RF - RS 2468 371 161 6.657 <.0001
##
## Results are averaged over the levels of: Treatment2
## Degrees-of-freedom method: containmentlibrary(Rmisc)
g12304_sum<-summarySE(biomin_all_filtered_long_g12304 , measurevar='Expression', groupvars=c('Origin', 'Treatment2'), na.rm=TRUE, conf.interval = 0.95)
g12304_sum## Origin Treatment2 N Expression sd se ci
## 1 RF Stable 44 5960.191 2423.309 365.3276 736.7533
## 2 RF Variable 44 4766.484 2433.696 366.8934 739.9111
## 3 RS Stable 48 2494.782 1512.223 218.2706 439.1038
## 4 RS Variable 48 2791.885 1392.250 200.9540 404.2674Figure
pd<- position_dodge(0.2)
g12304_fig<-ggplot(data=g12304_sum, aes(y=Expression, x=Treatment2, color=Origin),group = interaction(Origin))+
geom_point(data=biomin_all_filtered_long_g12304,aes(y=Expression, x=Treatment2, color=Origin), alpha=0.4, position = pd)+
geom_line(aes(group = interaction(Origin), stat="identity"),position=position_dodge(0.2))+
geom_point(size=3, stat="identity", position = pd)+
geom_errorbar(aes(ymin=Expression-se, ymax=Expression+se), stat="identity",width=0.2, position = pd)+
scale_fill_manual("Origin", values=c("RS"='#26519B', "RF"= "#FE180C"))+
scale_color_manual("Origin", values=c("RS"='#26519B', "RF"= "#FE180C"))+
scale_y_continuous(expression(Carbonic~anhydrase~("STPCA2-1")), limits=c(0,10000))+
ggtitle(~blue)+
theme_classic()+
theme(axis.text.x=element_text(vjust=0.5,size=12),#angling the labels on the x-axis
plot.title = element_text(margin = margin(t = 10, b = 10), hjust=0.5),#telling it where to position our plot title
panel.background= element_rect(fill=NA, color='black'),#this is making the black box around the graph
#strip.background = element_blank(),
#strip.text = element_blank(),
legend.title = element_text(vjust=0.5,size=12),
legend.position="none",
axis.text.y = element_text(vjust=0.5, size=12), #making the axis text larger
axis.title.x = element_blank(),#making the axis title larger
axis.title.y = element_text(size=12))#making the axis title larger ## Warning in geom_line(aes(group = interaction(Origin), stat = "identity"), :
## Ignoring unknown aesthetics: statg12304_fig## Warning: Removed 4 rows containing missing values (`geom_point()`).
SLC4 gamma- solute carrier family 4 member gamma
biomin_all_filtered_long_g15280 <- biomin_all_filtered_long %>%
filter(Gene == "Pocillopora_acuta_HIv2___RNAseq.g15280.t1")
biomin_all_filtered_long_g15280 ## # A tibble: 46 × 36
## Gene Dispersion AIC logLik meanExp X.Intercept. TreatmentVariable
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Pocillopora_a… 0.376 425. -206. 5.12 3.45 0.0890
## 2 Pocillopora_a… 0.376 425. -206. 5.12 3.45 0.0890
## 3 Pocillopora_a… 0.376 425. -206. 5.12 3.45 0.0890
## 4 Pocillopora_a… 0.376 425. -206. 5.12 3.45 0.0890
## 5 Pocillopora_a… 0.376 425. -206. 5.12 3.45 0.0890
## 6 Pocillopora_a… 0.376 425. -206. 5.12 3.45 0.0890
## 7 Pocillopora_a… 0.376 425. -206. 5.12 3.45 0.0890
## 8 Pocillopora_a… 0.376 425. -206. 5.12 3.45 0.0890
## 9 Pocillopora_a… 0.376 425. -206. 5.12 3.45 0.0890
## 10 Pocillopora_a… 0.376 425. -206. 5.12 3.45 0.0890
## # ℹ 36 more rows
## # ℹ 29 more variables: OriginFlat <dbl>, TreatmentVariable.OriginFlat <dbl>,
## # se_.Intercept. <dbl>, se_TreatmentVariable <dbl>, se_OriginFlat <dbl>,
## # se_TreatmentVariable.OriginFlat <dbl>, Chisq_Treatment <dbl>,
## # Chisq_Origin <dbl>, Chisq_Treatment.Origin <dbl>, Df_Treatment <int>,
## # Df_Origin <int>, Df_Treatment.Origin <int>, P_Treatment <dbl>,
## # P_Origin <dbl>, P_Treatment.Origin <dbl>, Treatment <fct>, …library(nlme)
g15280.lme <- lme(Expression~Origin*Treatment2, random = ~1|Colony, data=biomin_all_filtered_long_g15280 , na.action=na.exclude)
car::Anova(g15280.lme, type=3)## Analysis of Deviance Table (Type III tests)
##
## Response: Expression
## Chisq Df Pr(>Chisq)
## (Intercept) 82.1177 1 < 2.2e-16 ***
## Origin 9.0855 1 0.002576 **
## Treatment2 0.7652 1 0.381698
## Origin:Treatment2 0.9266 1 0.335743
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1tukey3<- emmeans(g15280.lme, list(pairwise ~ Origin), adjust = "tukey")## NOTE: Results may be misleading due to involvement in interactionstukey3## $`emmeans of Origin`
## Origin emmean SE df lower.CL upper.CL
## RF 52.0 4.36 19 42.9 61.2
## RS 32.1 4.17 19 23.4 40.9
##
## Results are averaged over the levels of: Treatment2
## Degrees-of-freedom method: containment
## Confidence level used: 0.95
##
## $`pairwise differences of Origin`
## 1 estimate SE df t.ratio p.value
## RF - RS 19.9 6.03 23 3.300 0.0031
##
## Results are averaged over the levels of: Treatment2
## Degrees-of-freedom method: containmentlibrary(Rmisc)
g15280_sum<-summarySE(biomin_all_filtered_long_g15280 , measurevar='Expression', groupvars=c('Origin', 'Treatment2'), na.rm=TRUE, conf.interval = 0.95)
g15280_sum## Origin Treatment2 N Expression sd se ci
## 1 RF Stable 11 55.84385 24.85863 7.495159 16.700254
## 2 RF Variable 11 48.22013 24.68527 7.442889 16.583790
## 3 RS Stable 12 30.12777 14.24152 4.111172 9.048629
## 4 RS Variable 12 34.11841 16.62673 4.799723 10.564118Figure
pd<- position_dodge(0.2)
g15280_fig<-ggplot(data=g15280_sum, aes(y=Expression, x=Treatment2, color=Origin),group = interaction(Origin))+
geom_point(data=biomin_all_filtered_long_g15280,aes(y=Expression, x=Treatment2, color=Origin), alpha=0.4, position = pd)+
geom_line(aes(group = interaction(Origin), stat="identity"),position=position_dodge(0.2))+
geom_point(size=3, stat="identity", position = pd)+
geom_errorbar(aes(ymin=Expression-se, ymax=Expression+se), stat="identity",width=0.2, position = pd)+
scale_fill_manual("Origin", values=c("RS"='#26519B', "RF"= "#FE180C"))+
scale_color_manual("Origin", values=c("RS"='#26519B', "RF"= "#FE180C"))+
scale_y_continuous(expression(Solute~carrier~family~4~member~gamma))+
ggtitle(~blue)+
theme_classic()+
theme(axis.text.x=element_text(vjust=0.5,size=12),#angling the labels on the x-axis
plot.title = element_text(margin = margin(t = 10, b = 10), hjust=0.5),#telling it where to position our plot title
panel.background= element_rect(fill=NA, color='black'),#this is making the black box around the graph
#strip.background = element_blank(),
#strip.text = element_blank(),
legend.title = element_text(vjust=0.5,size=12),
legend.position="none",
axis.text.y = element_text(vjust=0.5, size=12), #making the axis text larger
axis.title.x = element_blank(),#making the axis title larger
axis.title.y = element_text(size=12))#making the axis title larger ## Warning in geom_line(aes(group = interaction(Origin), stat = "identity"), :
## Ignoring unknown aesthetics: statg15280_fig
SLC4A7 - sodium bicarbonate cotransporter 3-like isoform X2
biomin_all_filtered_long_g7402 <- biomin_all_filtered_long %>%
filter(Gene == "Pocillopora_acuta_HIv2___RNAseq.g7402.t1")
biomin_all_filtered_long_g7402 ## # A tibble: 46 × 36
## Gene Dispersion AIC logLik meanExp X.Intercept. TreatmentVariable
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Pocillopora_a… 0.229 640. -314. 8.97 6.14 0.0129
## 2 Pocillopora_a… 0.229 640. -314. 8.97 6.14 0.0129
## 3 Pocillopora_a… 0.229 640. -314. 8.97 6.14 0.0129
## 4 Pocillopora_a… 0.229 640. -314. 8.97 6.14 0.0129
## 5 Pocillopora_a… 0.229 640. -314. 8.97 6.14 0.0129
## 6 Pocillopora_a… 0.229 640. -314. 8.97 6.14 0.0129
## 7 Pocillopora_a… 0.229 640. -314. 8.97 6.14 0.0129
## 8 Pocillopora_a… 0.229 640. -314. 8.97 6.14 0.0129
## 9 Pocillopora_a… 0.229 640. -314. 8.97 6.14 0.0129
## 10 Pocillopora_a… 0.229 640. -314. 8.97 6.14 0.0129
## # ℹ 36 more rows
## # ℹ 29 more variables: OriginFlat <dbl>, TreatmentVariable.OriginFlat <dbl>,
## # se_.Intercept. <dbl>, se_TreatmentVariable <dbl>, se_OriginFlat <dbl>,
## # se_TreatmentVariable.OriginFlat <dbl>, Chisq_Treatment <dbl>,
## # Chisq_Origin <dbl>, Chisq_Treatment.Origin <dbl>, Df_Treatment <int>,
## # Df_Origin <int>, Df_Treatment.Origin <int>, P_Treatment <dbl>,
## # P_Origin <dbl>, P_Treatment.Origin <dbl>, Treatment <fct>, …library(nlme)
g7402.lme <- lme(Expression~Origin*Treatment2, random = ~1|Colony, data=biomin_all_filtered_long_g7402 , na.action=na.exclude)
car::Anova(g7402.lme, type=3)## Analysis of Deviance Table (Type III tests)
##
## Response: Expression
## Chisq Df Pr(>Chisq)
## (Intercept) 87.9465 1 <2e-16 ***
## Origin 6.3555 1 0.0117 *
## Treatment2 0.8936 1 0.3445
## Origin:Treatment2 0.8273 1 0.3630
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1tukey3<- emmeans(g7402.lme, list(pairwise ~ Origin), adjust = "tukey")## NOTE: Results may be misleading due to involvement in interactionstukey3## $`emmeans of Origin`
## Origin emmean SE df lower.CL upper.CL
## RF 660 54.9 19 545 775
## RS 464 52.6 19 353 574
##
## Results are averaged over the levels of: Treatment2
## Degrees-of-freedom method: containment
## Confidence level used: 0.95
##
## $`pairwise differences of Origin`
## 1 estimate SE df t.ratio p.value
## RF - RS 197 75 23 2.621 0.0153
##
## Results are averaged over the levels of: Treatment2
## Degrees-of-freedom method: containmentlibrary(Rmisc)
g7402_sum<-summarySE(biomin_all_filtered_long_g7402 , measurevar='Expression', groupvars=c('Origin', 'Treatment2'), na.rm=TRUE, conf.interval = 0.95)
g7402_sum## Origin Treatment2 N Expression sd se ci
## 1 RF Stable 11 710.2028 264.1958 79.65802 177.4891
## 2 RF Variable 11 611.9007 312.5396 94.23425 209.9670
## 3 RS Stable 12 445.6090 228.0413 65.82984 144.8905
## 4 RS Variable 12 478.2523 190.4959 54.99144 121.0353Figure
pd<- position_dodge(0.2)
g7402_fig<-ggplot(data=g7402_sum, aes(y=Expression, x=Treatment2, color=Origin),group = interaction(Origin))+
geom_point(data=biomin_all_filtered_long_g7402,aes(y=Expression, x=Treatment2, color=Origin), alpha=0.4, position = pd)+
geom_line(aes(group = interaction(Origin), stat="identity"),position=position_dodge(0.2))+
geom_point(size=3, stat="identity", position = pd)+
geom_errorbar(aes(ymin=Expression-se, ymax=Expression+se), stat="identity",width=0.2, position = pd)+
scale_fill_manual("Origin", values=c("RS"='#26519B', "RF"= "#FE180C"))+
scale_color_manual("Origin", values=c("RS"='#26519B', "RF"= "#FE180C"))+
scale_y_continuous(expression(sodium~bicarbonate~cotransporter~3-like~isoform~X2))+
ggtitle(~blue)+
theme_classic()+
theme(axis.text.x=element_text(vjust=0.5,size=12),#angling the labels on the x-axis
plot.title = element_text(margin = margin(t = 10, b = 10), hjust=0.5),#telling it where to position our plot title
panel.background= element_rect(fill=NA, color='black'),#this is making the black box around the graph
#strip.background = element_blank(),
#strip.text = element_blank(),
legend.title = element_text(vjust=0.5,size=12),
legend.position="none",
axis.text.y = element_text(vjust=0.5, size=12), #making the axis text larger
axis.title.x = element_blank(),#making the axis title larger
axis.title.y = element_text(size=12))#making the axis title larger ## Warning in geom_line(aes(group = interaction(Origin), stat = "identity"), :
## Ignoring unknown aesthetics: statg7402_fig
Protein lingerer-like
biomin_all_filtered_long_g7908 <- biomin_all_filtered_long %>%
filter(Gene == "Pocillopora_acuta_HIv2___RNAseq.g7908.t1")
biomin_all_filtered_long_g7908 ## # A tibble: 46 × 36
## Gene Dispersion AIC logLik meanExp X.Intercept. TreatmentVariable
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Pocillopora_a… 0.138 426. -207. 5.95 3.98 -0.0743
## 2 Pocillopora_a… 0.138 426. -207. 5.95 3.98 -0.0743
## 3 Pocillopora_a… 0.138 426. -207. 5.95 3.98 -0.0743
## 4 Pocillopora_a… 0.138 426. -207. 5.95 3.98 -0.0743
## 5 Pocillopora_a… 0.138 426. -207. 5.95 3.98 -0.0743
## 6 Pocillopora_a… 0.138 426. -207. 5.95 3.98 -0.0743
## 7 Pocillopora_a… 0.138 426. -207. 5.95 3.98 -0.0743
## 8 Pocillopora_a… 0.138 426. -207. 5.95 3.98 -0.0743
## 9 Pocillopora_a… 0.138 426. -207. 5.95 3.98 -0.0743
## 10 Pocillopora_a… 0.138 426. -207. 5.95 3.98 -0.0743
## # ℹ 36 more rows
## # ℹ 29 more variables: OriginFlat <dbl>, TreatmentVariable.OriginFlat <dbl>,
## # se_.Intercept. <dbl>, se_TreatmentVariable <dbl>, se_OriginFlat <dbl>,
## # se_TreatmentVariable.OriginFlat <dbl>, Chisq_Treatment <dbl>,
## # Chisq_Origin <dbl>, Chisq_Treatment.Origin <dbl>, Df_Treatment <int>,
## # Df_Origin <int>, Df_Treatment.Origin <int>, P_Treatment <dbl>,
## # P_Origin <dbl>, P_Treatment.Origin <dbl>, Treatment <fct>, …library(nlme)
g7908.lme <- lme(Expression~Origin*Treatment2, random = ~1|Colony, data=biomin_all_filtered_long_g7908 , na.action=na.exclude)
car::Anova(g7908.lme, type=3)## Analysis of Deviance Table (Type III tests)
##
## Response: Expression
## Chisq Df Pr(>Chisq)
## (Intercept) 207.4212 1 < 2.2e-16 ***
## Origin 15.8613 1 6.816e-05 ***
## Treatment2 0.1605 1 0.6887
## Origin:Treatment2 0.0058 1 0.9392
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1tukey3<- emmeans(g7908.lme, list(pairwise ~ Origin), adjust = "tukey")## NOTE: Results may be misleading due to involvement in interactionstukey3## $`emmeans of Origin`
## Origin emmean SE df lower.CL upper.CL
## RF 81.2 4.07 19 72.7 89.7
## RS 49.9 3.89 19 41.8 58.0
##
## Results are averaged over the levels of: Treatment2
## Degrees-of-freedom method: containment
## Confidence level used: 0.95
##
## $`pairwise differences of Origin`
## 1 estimate SE df t.ratio p.value
## RF - RS 31.3 5.63 23 5.556 <.0001
##
## Results are averaged over the levels of: Treatment2
## Degrees-of-freedom method: containmentlibrary(Rmisc)
g7908_sum<-summarySE(biomin_all_filtered_long_g7908 , measurevar='Expression', groupvars=c('Origin', 'Treatment2'), na.rm=TRUE, conf.interval = 0.95)
g7908_sum## Origin Treatment2 N Expression sd se ci
## 1 RF Stable 11 82.80022 22.33628 6.734642 15.005717
## 2 RF Variable 11 79.54283 21.44921 6.467180 14.409774
## 3 RS Stable 12 51.10111 16.98594 4.903419 10.792352
## 4 RS Variable 12 48.70220 15.09645 4.357969 9.591824Figure
pd<- position_dodge(0.2)
g7908_fig<-ggplot(data=g7908_sum, aes(y=Expression, x=Treatment2, color=Origin),group = interaction(Origin))+
geom_point(data=biomin_all_filtered_long_g7908,aes(y=Expression, x=Treatment2, color=Origin), alpha=0.4, position = pd)+
geom_line(aes(group = interaction(Origin), stat="identity"),position=position_dodge(0.2))+
geom_point(size=3, stat="identity", position = pd)+
geom_errorbar(aes(ymin=Expression-se, ymax=Expression+se), stat="identity",width=0.2, position = pd)+
scale_fill_manual("Origin", values=c("RS"='#26519B', "RF"= "#FE180C"))+
scale_color_manual("Origin", values=c("RS"='#26519B', "RF"= "#FE180C"))+
scale_y_continuous(expression(protein~lingerer-like))+
ggtitle(~blue)+
theme_classic()+
theme(axis.text.x=element_text(vjust=0.5,size=12),#angling the labels on the x-axis
plot.title = element_text(margin = margin(t = 10, b = 10), hjust=0.5),#telling it where to position our plot title
panel.background= element_rect(fill=NA, color='black'),#this is making the black box around the graph
#strip.background = element_blank(),
#strip.text = element_blank(),
legend.title = element_text(vjust=0.5,size=12),
legend.position="none",
axis.text.y = element_text(vjust=0.5, size=12), #making the axis text larger
axis.title.x = element_blank(),#making the axis title larger
axis.title.y = element_text(size=12))#making the axis title larger ## Warning in geom_line(aes(group = interaction(Origin), stat = "identity"), :
## Ignoring unknown aesthetics: statg7908_fig
Signficant genes not in WGCNA modules
Vitellogenin
biomin_all_filtered_long_g13222 <- biomin_all_filtered_long %>%
filter(Gene == "Pocillopora_acuta_HIv2___TS.g13222.t1b")
biomin_all_filtered_long_g13222 ## # A tibble: 184 × 36
## Gene Dispersion AIC logLik meanExp X.Intercept. TreatmentVariable
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Pocillopora_a… 0.0506 621. -305. 9.50 6.53 -0.123
## 2 Pocillopora_a… 0.0506 621. -305. 9.50 6.53 -0.123
## 3 Pocillopora_a… 0.0506 621. -305. 9.50 6.53 -0.123
## 4 Pocillopora_a… 0.0506 621. -305. 9.50 6.53 -0.123
## 5 Pocillopora_a… 0.0506 621. -305. 9.50 6.53 -0.123
## 6 Pocillopora_a… 0.0506 621. -305. 9.50 6.53 -0.123
## 7 Pocillopora_a… 0.0506 621. -305. 9.50 6.53 -0.123
## 8 Pocillopora_a… 0.0506 621. -305. 9.50 6.53 -0.123
## 9 Pocillopora_a… 0.0506 621. -305. 9.50 6.53 -0.123
## 10 Pocillopora_a… 0.0506 621. -305. 9.50 6.53 -0.123
## # ℹ 174 more rows
## # ℹ 29 more variables: OriginFlat <dbl>, TreatmentVariable.OriginFlat <dbl>,
## # se_.Intercept. <dbl>, se_TreatmentVariable <dbl>, se_OriginFlat <dbl>,
## # se_TreatmentVariable.OriginFlat <dbl>, Chisq_Treatment <dbl>,
## # Chisq_Origin <dbl>, Chisq_Treatment.Origin <dbl>, Df_Treatment <int>,
## # Df_Origin <int>, Df_Treatment.Origin <int>, P_Treatment <dbl>,
## # P_Origin <dbl>, P_Treatment.Origin <dbl>, Treatment <fct>, …library(nlme)
g13222.lme <- lme(Expression~Origin*Treatment2, random = ~1|Colony, data=biomin_all_filtered_long_g13222 , na.action=na.exclude)
car::Anova(g13222.lme, type=3)## Analysis of Deviance Table (Type III tests)
##
## Response: Expression
## Chisq Df Pr(>Chisq)
## (Intercept) 1126.7286 1 < 2.2e-16 ***
## Origin 41.2035 1 1.372e-10 ***
## Treatment2 2.9349 1 0.08669 .
## Origin:Treatment2 0.5477 1 0.45926
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1tukey3<- emmeans(g13222.lme, list(pairwise ~ Origin), adjust = "tukey")## NOTE: Results may be misleading due to involvement in interactionstukey3## $`emmeans of Origin`
## Origin emmean SE df lower.CL upper.CL
## RF 832 22.8 19 784 879
## RS 636 22.2 19 590 683
##
## Results are averaged over the levels of: Treatment2
## Degrees-of-freedom method: containment
## Confidence level used: 0.95
##
## $`pairwise differences of Origin`
## 1 estimate SE df t.ratio p.value
## RF - RS 195 24.2 161 8.053 <.0001
##
## Results are averaged over the levels of: Treatment2
## Degrees-of-freedom method: containmentlibrary(Rmisc)
g13222_sum<-summarySE(biomin_all_filtered_long_g13222 , measurevar='Expression', groupvars=c('Origin', 'Treatment2'), na.rm=TRUE, conf.interval = 0.95)
g13222_sum## Origin Treatment2 N Expression sd se ci
## 1 RF Stable 44 859.5886 158.58611 23.90776 48.21459
## 2 RF Variable 44 821.8027 114.92648 17.32582 34.94084
## 3 RS Stable 48 660.1489 146.62900 21.16407 42.57662
## 4 RS Variable 48 599.7645 70.60441 10.19087 20.50138Figure
pd<- position_dodge(0.2)
g13222_fig<-ggplot(data=g13222_sum, aes(y=Expression, x=Treatment2, color=Origin),group = interaction(Origin))+
geom_point(data=biomin_all_filtered_long_g13222,aes(y=Expression, x=Treatment2, color=Origin), alpha=0.4, position = pd)+
geom_line(aes(group = interaction(Origin), stat="identity"),position=position_dodge(0.2))+
geom_point(size=3, stat="identity", position = pd)+
geom_errorbar(aes(ymin=Expression-se, ymax=Expression+se), stat="identity",width=0.2, position = pd)+
scale_fill_manual("Origin", values=c("RS"='#26519B', "RF"= "#FE180C"))+
scale_color_manual("Origin", values=c("RS"='#26519B', "RF"= "#FE180C"))+
scale_y_continuous(expression(Vitellogenin))+
ggtitle(~blue)+
theme_classic()+
theme(axis.text.x=element_text(vjust=0.5,size=12),#angling the labels on the x-axis
plot.title = element_text(margin = margin(t = 10, b = 10), hjust=0.5),#telling it where to position our plot title
panel.background= element_rect(fill=NA, color='black'),#this is making the black box around the graph
#strip.background = element_blank(),
#strip.text = element_blank(),
legend.title = element_text(vjust=0.5,size=12),
legend.position="none",
axis.text.y = element_text(vjust=0.5, size=12), #making the axis text larger
axis.title.x = element_blank(),#making the axis title larger
axis.title.y = element_text(size=12))#making the axis title larger ## Warning in geom_line(aes(group = interaction(Origin), stat = "identity"), :
## Ignoring unknown aesthetics: statg13222_fig
uncharacterized protein LOC111323869
biomin_all_filtered_long_g21232 <- biomin_all_filtered_long %>%
filter(Gene == "Pocillopora_acuta_HIv2___RNAseq.g21232.t1")
biomin_all_filtered_long_g21232 ## # A tibble: 46 × 36
## Gene Dispersion AIC logLik meanExp X.Intercept. TreatmentVariable
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Pocillopora_a… 0.0605 353. -171. 5.22 3.54 -0.0600
## 2 Pocillopora_a… 0.0605 353. -171. 5.22 3.54 -0.0600
## 3 Pocillopora_a… 0.0605 353. -171. 5.22 3.54 -0.0600
## 4 Pocillopora_a… 0.0605 353. -171. 5.22 3.54 -0.0600
## 5 Pocillopora_a… 0.0605 353. -171. 5.22 3.54 -0.0600
## 6 Pocillopora_a… 0.0605 353. -171. 5.22 3.54 -0.0600
## 7 Pocillopora_a… 0.0605 353. -171. 5.22 3.54 -0.0600
## 8 Pocillopora_a… 0.0605 353. -171. 5.22 3.54 -0.0600
## 9 Pocillopora_a… 0.0605 353. -171. 5.22 3.54 -0.0600
## 10 Pocillopora_a… 0.0605 353. -171. 5.22 3.54 -0.0600
## # ℹ 36 more rows
## # ℹ 29 more variables: OriginFlat <dbl>, TreatmentVariable.OriginFlat <dbl>,
## # se_.Intercept. <dbl>, se_TreatmentVariable <dbl>, se_OriginFlat <dbl>,
## # se_TreatmentVariable.OriginFlat <dbl>, Chisq_Treatment <dbl>,
## # Chisq_Origin <dbl>, Chisq_Treatment.Origin <dbl>, Df_Treatment <int>,
## # Df_Origin <int>, Df_Treatment.Origin <int>, P_Treatment <dbl>,
## # P_Origin <dbl>, P_Treatment.Origin <dbl>, Treatment <fct>, …library(nlme)
g21232.lme <- lme(Expression~Origin*Treatment2, random = ~1|Colony, data=biomin_all_filtered_long_g21232 , na.action=na.exclude)
car::Anova(g21232.lme, type=3)## Analysis of Deviance Table (Type III tests)
##
## Response: Expression
## Chisq Df Pr(>Chisq)
## (Intercept) 209.9697 1 < 2e-16 ***
## Origin 4.7421 1 0.02943 *
## Treatment2 0.1307 1 0.71769
## Origin:Treatment2 0.0003 1 0.98700
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1tukey3<- emmeans(g21232.lme, list(pairwise ~ Origin), adjust = "tukey")## NOTE: Results may be misleading due to involvement in interactionstukey3## $`emmeans of Origin`
## Origin emmean SE df lower.CL upper.CL
## RF 42.0 2.4 19 36.9 47.0
## RS 33.7 2.3 19 28.8 38.5
##
## Results are averaged over the levels of: Treatment2
## Degrees-of-freedom method: containment
## Confidence level used: 0.95
##
## $`pairwise differences of Origin`
## 1 estimate SE df t.ratio p.value
## RF - RS 8.29 3.02 23 2.747 0.0115
##
## Results are averaged over the levels of: Treatment2
## Degrees-of-freedom method: containmentlibrary(Rmisc)
g21232_sum<-summarySE(biomin_all_filtered_long_g21232 , measurevar='Expression', groupvars=c('Origin', 'Treatment2'), na.rm=TRUE, conf.interval = 0.95)
g21232_sum## Origin Treatment2 N Expression sd se ci
## 1 RF Stable 11 42.29276 10.185331 3.070993 6.842599
## 2 RF Variable 11 41.06536 11.600467 3.497672 7.793300
## 3 RS Stable 12 33.84473 7.075632 2.042559 4.495642
## 4 RS Variable 12 32.69394 10.921097 3.152649 6.938934Figure
pd<- position_dodge(0.2)
g21232_fig<-ggplot(data=g21232_sum, aes(y=Expression, x=Treatment2, color=Origin),group = interaction(Origin))+
geom_point(data=biomin_all_filtered_long_g21232,aes(y=Expression, x=Treatment2, color=Origin), alpha=0.4, position = pd)+
geom_line(aes(group = interaction(Origin), stat="identity"),position=position_dodge(0.2))+
geom_point(size=3, stat="identity", position = pd)+
geom_errorbar(aes(ymin=Expression-se, ymax=Expression+se), stat="identity",width=0.2, position = pd)+
scale_fill_manual("Origin", values=c("RS"='#26519B', "RF"= "#FE180C"))+
scale_color_manual("Origin", values=c("RS"='#26519B', "RF"= "#FE180C"))+
scale_y_continuous(expression(uncharacterized~protein~LOC111323869))+
ggtitle(~blue)+
theme_classic()+
theme(axis.text.x=element_text(vjust=0.5,size=12),#angling the labels on the x-axis
plot.title = element_text(margin = margin(t = 10, b = 10), hjust=0.5),#telling it where to position our plot title
panel.background= element_rect(fill=NA, color='black'),#this is making the black box around the graph
#strip.background = element_blank(),
#strip.text = element_blank(),
legend.title = element_text(vjust=0.5,size=12),
legend.position="none",
axis.text.y = element_text(vjust=0.5, size=12), #making the axis text larger
axis.title.x = element_blank(),#making the axis title larger
axis.title.y = element_text(size=12))#making the axis title larger ## Warning in geom_line(aes(group = interaction(Origin), stat = "identity"), :
## Ignoring unknown aesthetics: statg21232_fig
uncharacterized protein LOC111345150
biomin_all_filtered_long_g20587 <- biomin_all_filtered_long %>%
filter(Gene == "Pocillopora_acuta_HIv2___RNAseq.g20587.t2")
biomin_all_filtered_long_g20587 ## # A tibble: 46 × 36
## Gene Dispersion AIC logLik meanExp X.Intercept. TreatmentVariable
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Pocillopora_a… 4.85 271. -129. 1.93 2.59 0.106
## 2 Pocillopora_a… 4.85 271. -129. 1.93 2.59 0.106
## 3 Pocillopora_a… 4.85 271. -129. 1.93 2.59 0.106
## 4 Pocillopora_a… 4.85 271. -129. 1.93 2.59 0.106
## 5 Pocillopora_a… 4.85 271. -129. 1.93 2.59 0.106
## 6 Pocillopora_a… 4.85 271. -129. 1.93 2.59 0.106
## 7 Pocillopora_a… 4.85 271. -129. 1.93 2.59 0.106
## 8 Pocillopora_a… 4.85 271. -129. 1.93 2.59 0.106
## 9 Pocillopora_a… 4.85 271. -129. 1.93 2.59 0.106
## 10 Pocillopora_a… 4.85 271. -129. 1.93 2.59 0.106
## # ℹ 36 more rows
## # ℹ 29 more variables: OriginFlat <dbl>, TreatmentVariable.OriginFlat <dbl>,
## # se_.Intercept. <dbl>, se_TreatmentVariable <dbl>, se_OriginFlat <dbl>,
## # se_TreatmentVariable.OriginFlat <dbl>, Chisq_Treatment <dbl>,
## # Chisq_Origin <dbl>, Chisq_Treatment.Origin <dbl>, Df_Treatment <int>,
## # Df_Origin <int>, Df_Treatment.Origin <int>, P_Treatment <dbl>,
## # P_Origin <dbl>, P_Treatment.Origin <dbl>, Treatment <fct>, …library(nlme)
g20587.lme <- lme(Expression~Origin*Treatment2, random = ~1|Colony, data=biomin_all_filtered_long_g20587 , na.action=na.exclude)
car::Anova(g20587.lme, type=3)## Analysis of Deviance Table (Type III tests)
##
## Response: Expression
## Chisq Df Pr(>Chisq)
## (Intercept) 0.1659 1 0.683767
## Origin 7.2048 1 0.007271 **
## Treatment2 0.1684 1 0.681502
## Origin:Treatment2 0.0040 1 0.949524
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1tukey3<- emmeans(g20587.lme, list(pairwise ~ Origin), adjust = "tukey")## NOTE: Results may be misleading due to involvement in interactionstukey3## $`emmeans of Origin`
## Origin emmean SE df lower.CL upper.CL
## RF 2.06 2.78 19 -3.77 7.88
## RS 13.31 2.68 19 7.70 18.93
##
## Results are averaged over the levels of: Treatment2
## Degrees-of-freedom method: containment
## Confidence level used: 0.95
##
## $`pairwise differences of Origin`
## 1 estimate SE df t.ratio p.value
## RF - RS -11.3 3.35 23 -3.355 0.0027
##
## Results are averaged over the levels of: Treatment2
## Degrees-of-freedom method: containmentlibrary(Rmisc)
g20587_sum<-summarySE(biomin_all_filtered_long_g20587 , measurevar='Expression', groupvars=c('Origin', 'Treatment2'), na.rm=TRUE, conf.interval = 0.95)
g20587_sum## Origin Treatment2 N Expression sd se ci
## 1 RF Stable 11 1.069422 2.379314 0.7173901 1.598445
## 2 RF Variable 11 2.503911 3.677310 1.1087507 2.470451
## 3 RS Stable 12 12.756815 15.369960 4.4369253 9.765607
## 4 RS Variable 12 14.497631 15.007069 4.3321675 9.535036Figure
pd<- position_dodge(0.2)
g20587_fig<-ggplot(data=g20587_sum, aes(y=Expression, x=Treatment2, color=Origin),group = interaction(Origin))+
geom_point(data=biomin_all_filtered_long_g20587,aes(y=Expression, x=Treatment2, color=Origin), alpha=0.4, position = pd)+
geom_line(aes(group = interaction(Origin), stat="identity"),position=position_dodge(0.2))+
geom_point(size=3, stat="identity", position = pd)+
geom_errorbar(aes(ymin=Expression-se, ymax=Expression+se), stat="identity",width=0.2, position = pd)+
scale_fill_manual("Origin", values=c("RS"='#26519B', "RF"= "#FE180C"))+
scale_color_manual("Origin", values=c("RS"='#26519B', "RF"= "#FE180C"))+
scale_y_continuous(expression(uncharacterized~protein~LOC111345150), limits=c(0,30))+
ggtitle(~blue)+
theme_classic()+
theme(axis.text.x=element_text(vjust=0.5,size=12),#angling the labels on the x-axis
plot.title = element_text(margin = margin(t = 10, b = 10), hjust=0.5),#telling it where to position our plot title
panel.background= element_rect(fill=NA, color='black'),#this is making the black box around the graph
#strip.background = element_blank(),
#strip.text = element_blank(),
legend.title = element_text(vjust=0.5,size=12),
legend.position="none",
axis.text.y = element_text(vjust=0.5, size=12), #making the axis text larger
axis.title.x = element_blank(),#making the axis title larger
axis.title.y = element_text(size=12))#making the axis title larger ## Warning in geom_line(aes(group = interaction(Origin), stat = "identity"), :
## Ignoring unknown aesthetics: statg20587_fig## Warning: Removed 4 rows containing missing values (`geom_point()`).
Late embryogenesis protein
biomin_all_filtered_long_g16715 <- biomin_all_filtered_long %>%
filter(Gene == "Pocillopora_acuta_HIv2___RNAseq.g16715.t1")
biomin_all_filtered_long_g16715 ## # A tibble: 46 × 36
## Gene Dispersion AIC logLik meanExp X.Intercept. TreatmentVariable
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Pocillopora_a… 0.154 847. -417. 12.6 8.39 0.255
## 2 Pocillopora_a… 0.154 847. -417. 12.6 8.39 0.255
## 3 Pocillopora_a… 0.154 847. -417. 12.6 8.39 0.255
## 4 Pocillopora_a… 0.154 847. -417. 12.6 8.39 0.255
## 5 Pocillopora_a… 0.154 847. -417. 12.6 8.39 0.255
## 6 Pocillopora_a… 0.154 847. -417. 12.6 8.39 0.255
## 7 Pocillopora_a… 0.154 847. -417. 12.6 8.39 0.255
## 8 Pocillopora_a… 0.154 847. -417. 12.6 8.39 0.255
## 9 Pocillopora_a… 0.154 847. -417. 12.6 8.39 0.255
## 10 Pocillopora_a… 0.154 847. -417. 12.6 8.39 0.255
## # ℹ 36 more rows
## # ℹ 29 more variables: OriginFlat <dbl>, TreatmentVariable.OriginFlat <dbl>,
## # se_.Intercept. <dbl>, se_TreatmentVariable <dbl>, se_OriginFlat <dbl>,
## # se_TreatmentVariable.OriginFlat <dbl>, Chisq_Treatment <dbl>,
## # Chisq_Origin <dbl>, Chisq_Treatment.Origin <dbl>, Df_Treatment <int>,
## # Df_Origin <int>, Df_Treatment.Origin <int>, P_Treatment <dbl>,
## # P_Origin <dbl>, P_Treatment.Origin <dbl>, Treatment <fct>, …library(nlme)
g16715.lme <- lme(Expression~Origin*Treatment2, random = ~1|Colony, data=biomin_all_filtered_long_g16715 , na.action=na.exclude)
car::Anova(g16715.lme, type=3)## Analysis of Deviance Table (Type III tests)
##
## Response: Expression
## Chisq Df Pr(>Chisq)
## (Intercept) 265.1538 1 < 2.2e-16 ***
## Origin 32.7124 1 1.069e-08 ***
## Treatment2 0.3604 1 0.5483
## Origin:Treatment2 1.3445 1 0.2462
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1tukey3<- emmeans(g16715.lme, list(pairwise ~ Origin), adjust = "tukey")## NOTE: Results may be misleading due to involvement in interactionstukey3## $`emmeans of Origin`
## Origin emmean SE df lower.CL upper.CL
## RF 8304 381 19 7507 9101
## RS 4973 365 19 4209 5737
##
## Results are averaged over the levels of: Treatment2
## Degrees-of-freedom method: containment
## Confidence level used: 0.95
##
## $`pairwise differences of Origin`
## 1 estimate SE df t.ratio p.value
## RF - RS 3331 505 23 6.602 <.0001
##
## Results are averaged over the levels of: Treatment2
## Degrees-of-freedom method: containmentlibrary(Rmisc)
g16715_sum<-summarySE(biomin_all_filtered_long_g16715 , measurevar='Expression', groupvars=c('Origin', 'Treatment2'), na.rm=TRUE, conf.interval = 0.95)
g16715_sum## Origin Treatment2 N Expression sd se ci
## 1 RF Stable 11 8076.757 1731.637 522.1082 1163.3295
## 2 RF Variable 11 8462.423 1729.677 521.5172 1162.0128
## 3 RS Stable 12 4230.613 1879.259 542.4954 1194.0243
## 4 RS Variable 12 5647.539 1237.158 357.1369 786.0529Figure
pd<- position_dodge(0.2)
g16715_fig<-ggplot(data=g16715_sum, aes(y=Expression, x=Treatment2, color=Origin),group = interaction(Origin))+
geom_point(data=biomin_all_filtered_long_g16715,aes(y=Expression, x=Treatment2, color=Origin), alpha=0.4, position = pd)+
geom_line(aes(group = interaction(Origin), stat="identity"),position=position_dodge(0.2))+
geom_point(size=3, stat="identity", position = pd)+
geom_errorbar(aes(ymin=Expression-se, ymax=Expression+se), stat="identity",width=0.2, position = pd)+
scale_fill_manual("Origin", values=c("RS"='#26519B', "RF"= "#FE180C"))+
scale_color_manual("Origin", values=c("RS"='#26519B', "RF"= "#FE180C"))+
scale_y_continuous(expression(Late~embryogenesis~protein))+
ggtitle(~blue)+
theme_classic()+
theme(axis.text.x=element_text(vjust=0.5,size=12),#angling the labels on the x-axis
plot.title = element_text(margin = margin(t = 10, b = 10), hjust=0.5),#telling it where to position our plot title
panel.background= element_rect(fill=NA, color='black'),#this is making the black box around the graph
#strip.background = element_blank(),
#strip.text = element_blank(),
legend.title = element_text(vjust=0.5,size=12),
legend.position="none",
axis.text.y = element_text(vjust=0.5, size=12), #making the axis text larger
axis.title.x = element_blank(),#making the axis title larger
axis.title.y = element_text(size=12))#making the axis title larger ## Warning in geom_line(aes(group = interaction(Origin), stat = "identity"), :
## Ignoring unknown aesthetics: statg16715_fig
Comparison of all sig. genes
biomin_compare_figs<-cowplot::plot_grid(g13824_fig,g25351_fig,g10093_fig,g5013_fig, g12304_fig,g7908_fig,g12304_fig,g13222_fig,g16715_fig,g20587_fig,g21232_fig,g15280_fig,g7402_fig, nrow=4)## Warning: Removed 1 rows containing missing values (`geom_point()`).
## Removed 1 rows containing missing values (`geom_point()`).## Warning: Removed 4 rows containing missing values (`geom_point()`).
## Removed 4 rows containing missing values (`geom_point()`).
## Removed 4 rows containing missing values (`geom_point()`).biomin_compare_figs
ggsave(filename="../../output/WGCNA/GO_analysis/biomin_compare_figs.png", plot=biomin_compare_figs, dpi=300, height=10, units="in", limitsize=FALSE)## Saving 7 x 10 in image